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Systematic Analysis of Thermal Turbojet 
Propulsion 
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SUMMARY 


The propulsion unit is divisible into six sections: free inflow, 
intake, compressor, heat chamber, turbine, and exit nozzle. It 
is treated uniformly by the quintuple system of dynamic and 
thermodynamic flow equations with their power or heat inputs, 
the conditions of change of state (cycle), the continuity equa- 
tions, and the thrust equation. These equations also comprise 
in each section the retarding forces and the heat creation of 
turbulent flow and the loss of mechanical power and of heat 
leakage. They form successive sets of equations, and their 
solution is derived in terms of the required data of altitude, 
flying speed, thrust, pressures in front and at the exit of the 
compressor, permissible entrance temperature at the turbine, 
and permissible ratio of power of compressor and turbine to heat 
input. 

This procedure, to illustrate the basic ideas, is, at first, worked 
out for a perfect unit (that is, without losses), and is afterwards 
completed for the percentual losses of heat and work in each 
section. 

The necessary heat input, the propulsion effectiveness, the 
distribution of temperatures, pressures, velocities and cross 
sections, the best cycle conditions, and the most efficient pressure- 
velocity relations at the exit are then determined by those 
data. 


INTRODUCTION 


, ee THERMAL JET PROPULSION UNIT treated in this 
paper is given by a schematic representation of 
Fig. 1. It shows the five sections of intake, com- 
pressor, heat chamber with flame holders, turbine, 
and nozzle. The unit is supposed to move at an 
altitude (a,) defined by the density (p) and the 
pressure (fo) with the constant, not necessarily hori- 
zontal velocity (Vo), and is required to produce a thrust 
(Th.) sufficient to overcome the total drag of the 
airplane at a steady, horizontal, climbing, or gliding 
velocity (Vo). 


Presented at the Aerodynamics Session, Fourteenth Annual 
Meeting, I.A.S., New York, January 29-31, 1946. Revised and 
received December 16, 1946. 

* Professor of Aerodynamics and Aero Structures. 


The theoretical treatment of this paper was started 
in 1945 and presented at the Annual Meeting of the 
Institute of the Aeronautical Sciences in January, 
1946. Other work meanwhile has delayed the pub- 
lication, and this delay has been used to make some 
slight changes and additions. 

In view of the great number of valuable literary 
contributions in this field,t it may be appropriate to 
indicate the different approach and analysis of this 
paper. 

It is more elementary insofar as it starts from the 
Bernoulli equation of steady compressible flow, com- 
bined with the transformation first given by Zeuner, 
which expresses the influence of heat input and of dis- 
sipation of mechanic energy on the equation of flow. 

On the other hand, it attempts a more general analy- 
sis by incorporating in the integration of the flow 
equations the friction and leakage losses, the relations 
between change of state (cycle) and the changes of 
cross sections and the excess of turbine over compressor 
power. A further difference appears in the evaluation 
of the thrust of the jet in terms of pressure ratio, 
velocity ratio, and temperature ratio at the exit and in a 
discussion of a pertinent maximum condition. 

Compressor and turbine are represented by a force 


field, helical in each stage, which represents the actually 


+ See Bibliography at the end of this paper. 
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finite number of pressure-discontinuity surfaces of the 
blades by a continuous distribution of helical surfaces 
of infinitesimal discontinuity. This concept has the 
advantage that all duct sections can be treated by one 
and the same general flow equation and by its general 
integral and that they differ only in the values of the 
heat and mechanical work input or loss and the change 
of state characteristic. It can be proved that the 
actual compressor or turbine with its rotating blades 
and countervanes can be designed so as to accomplish 
the same effect as the force field.* 

The solution proceeds in successive steps from one 
section to the next in such a way that the requirements 
of compression ratio of the compressor, the temperature 
requirement in the heat chamber, and the thrust re- 
quirement determine the total heat input and the 
compressor and turbine power. 

The representation uses as much as possible di- 
mensionless variables and the energy and power values, 
including heat in mechanical units (ft.lb. per lb.) and 
(ft.lb. per lb.sec.) which unifies and simplifies the 
equations but can easily be changed into the language 
of B.t.u. or kg.cal. if so desired. 

An assumption of the analysis is that the approach 
of the relative flow velocity toward the velocity of sound 
does not give rise to noticeable shock waves in the 
inside of the duct. The velocity of sound being pro- 
portional to the square root of the temperature and the 
temperature ratio being between two and four every- 
where, except at the intake, gives a velocity of sound 
between 987 and 1,520 m.p.h. so that the approach 
to a local Mach Number 1 is not imminent if increases 
of local velocity in corners and bends are avoided. 
Only the intake temperatures will be smaller unless 
the intake is also heated by heat conduction from the 
compressor. In order to take into account the case of 
turbulent friction in the duct section in a not too 
complicated manner, the cross sections are assumed to 
be annular in such a measure as to enable a one-di- 
mensional treatment with an average velocity (super- 
posed by turbulence fluctuation), an average tempera- 
ture, and an average pressure in every cross section. 


SYMBOLS 


= pressure, density, absolute tem- 
perature, velocity relative to 
the propulsion duct, and area 
of a cross section, respectively 
= the dimensionless pressure dens- 
V ity, temperature, and velocity 
3 variables referred to the un- 
Vo disturbed atmosphere (Vy 
velocity of duct relative to 

atmosphere) 


b,p, T, V, A 


* The author has elaborated this streamline design method in 
a (hectographed) course on compressors and turbines at the 
Polytechnic Institute of Brooklyn in the 1945-1946 term. It 
was shown in this course how the transition from the continuous 
force and streamline field to a discontinuous pressure and 
streamline field can be accomplished. 


h — Vo?/2g 

H = po/"/200 Vo? = static over dynamic pressure 

Hh = po/pog = 

by; = lengths of the duct sections 

ds = element of length of a streamline 

ds’ =ds/l 

= hydraulic 

section 
relative loss of head per unit of 
length 

dQ/ds = heat input (if negative loss of 
heat) or heat created by dissi- 
pation of energy, per unit of 
length of travel and per unit 
of weight of a fluid element, 
in mechanical units (ft.lb per 
ft.Ib.) 

0= f;*d0 = heat input in mechanical units 
per unit of weight of fluid for 
the total length of travel along 
a duct section 

dW/ds = work input (if negative output) 
per unit of weight of fluid and 
per unit of length of travel (of 
a fluid element) 

W = Py dw = work input per unit of weight of 
fluid for the total length of 
travel (along a duct section) 


diameter along a 


Th. = thrust of jet propulsion unit 

Th.* = Th./poVoA,. = thrust coefficient 

m’ = pVA = mass passing duct per unit of 
time 

Owm'g = total external heat-power input 
per second ft.lb. per sec.) 

n = Th. Vo/Q.m'g = propulsive efficiency 

Y = ratio of specific heats of the gas 


in its actual state 


n = polytropic exponent 
a = n/(n — 1) 

Cu = p/(u — 

oN = A/(u — d) 

= ¢1/4D, 

a = 

w* = Wh- 

Q* = Qh-! 

b = a(W* — 40*)-1 fold X 

nS 

€ W./ 

= coefficient indicating the _per- 
centual loss of heat in a section, 
with the subscript pertinent 
to it 

Be = b(1 + 

Ba = crrba(1 — ba) 

= bi(1 — 

Ye = — 

Yh = (1 — 6, + Br)/Hp 

Yt = (1 — 6 — B:)/Hp 

Yn = 6n(1 + cynba)/Hu 

a = (1 — a)/Hp, with a subscript 
indicating the section 

Subscripts 

0 = undisturbed atmosphere at flying 
altitude 

t = duct section of intake and in 
general inlet of any section 

c = duct section of compressor 
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duct section of heat chamber 

duct section of turbine 

duct section of nozzle 

exit of a section and particularly 
of nozzle 


THE CONTINUITY EQUATION 


The continuity equation can be written in terms of 
the density p, the velocity V averaged over a cross 
section, and the cross section area A as functions of 
the distance s from the entrance of a duct-section. 


pVA = const. (1) 


Introducing the ratios of p and V relative to the values 
p, Vo in the undisturbed (ambient) atmosphere, viz., 


p/m =a, 


one can express Eq. (1) by 
ogA = Ay (la) 


where Ao is the cross section of the relative inflow far 
infront, where o = 1 (p = pp) andg = 1(V = Vp). 

In case the equation of state and the change of state 
is given, it will often be of advantage to express o in 
this Eq. (la) either by the temperature ratios 7)/T = + 
or by the pressure ratio p/p) = g 

Furthermore, if the flow is rotating so that the axial 
velocities are Visia = V cos 6 (@ pitch angle), Eq. 
(la) changes to og cos A = Ap. 


THE THRUST AND EFFICIENCY EQUATIONS 


Thrust and body drag of the propulsion unit are 
linked to each other by their influence on the surround- 
ing atmosphere. 

To express the performance of jet propulsion, one 
needs the separation of these two forces. This may 
be done in the following way: 

Consider a cylindrical fluid volume surrounding the 
propulsion unit and axially coinciding with the axis 
of the unit. The front section of this surrounding 
gas-fluid cylinder may be so far in front that on it 
one can assume undisturbed pressure pp and undis- 
turbed relative velocity Vo. Its rear cross section may 
coincide with the exit cross section of the unit but 
may radially continue to the cylindrical (generator) 
surface of the fluid cylinder, so far outside that only 
undisturbed pressures and velocities on this cylinder 
surface are encountered. 

However, the radially outside part of the rear cross 
section of said fluid cylinder, surrounding the duct unit, 
will have pressures and velocities disturbed by the 
outside drag of the duct body. This drag will be the 
drag of the outside of this body, while the drag of the 
entrance and the exit edges of it will have no effect, 
since they are streamlines by the incoming and outgoing 
jet. On the other hand, this consideration will also 
cover the eventual shock waves starting from the out- 
side or the front edges of the duct, since they only 
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contribute to the pressure and velocity distribution of 
the outside and the outside fluid cross section at the rear. 
Observing now the time rate of momentum of flow 
relative to the momentum of the undisturbed flow and 
the pressure at the exit section relative to the un- 
disturbed pressure of the atmosphere, it is possible to 
evaluate the thrust (Th.) of the jet propulsion duct, 
while for the evaluation of the body drag the observa- 
tion of momentum and pressure difference on the 
outside rear cross section of the surrounding fluid 
cylinder would be necessary. The thrust of the jet 
proper (subscript e for exit) is therefore given by: 


Th. = — Vo) + — fo) (2) 


or in the dimensionless variables defined above, by a 
thrust coefficient 


Th. Po 
Th.* = = ¢,¢,(¢ — 1 — 1 2 
— = 1) + 1) 


By observing that by the equation of state of the gas 
(b/po = RgpT/RgpoT. or g = or), o = and 
with the abbreviation H = po/'/2poVo?,. the above 
equation is changed to: 


Th.* = gre — 1) + — 1) (2b) 


Dividing the work per unit of time of the thrust, Eq. 
(2a), by the total heat input Q per unit of time and per 
unit of weight of the fluid times the weight of flow per 
unit of time, one arrives at an expression for the value 
of the propulsion efficiency (7) 


n aut [2(¢ )+ ge") 
(3) 


where h = Vo*/2g is the velocity head and Q ft.lb. per 
Ib. the heat input per unit of fluid weight in mechanical 
units. 

Eqs. (2), (2b), and (3) leave the question open as to 
how much velocity V,(¢,) and how much pressure 
b.(q-) should contribute to attain the maximum of 
thrust. The answer will depend on the obtainable 
pressure excess p, — and on the cross-section dis- 
tribution along the duct, and it will be given in one 
of the following sections. 

Another remark on the nature of the thrust may be 
clarifying. It was necessary to determine the thrust 
value by the indirect, somewhat abstract method of 
using a cylindrical air volume surrounding the pro- 
pulsion duct. The direct (concrete) method would 
derive the thrust as the axial component of the re- 
sultant force of all inside pressures and frictions on 
walls, blades, etc., of the duct unit. It is evident that 
such an evaluation, though correct in principle, is prac- 
tically impossible, This justifies the application of the 
somewhat indirect method used above of the excess of 
the leaving over the entering momentum and pressure 
components. 
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The separate determination of net thrust and ex- 
ternal drag, as explained above, is based on the premise 
that the outside surface of the duct does not contribute 
to the thrust. The author owes to J. V. Foa the 
interesting remark that a small contribution to the 
thrust may be produced if the inflow is partly spilling 
over the rounded intake edge and producing a suction, 
in case the intake is built such that a stagnation line 
appears on the inside near the edge. However, nothing 
seems to be gained because of the corresponding loss 
of inside mass flow, and at any rate not for compressor 
driven units. 


Tue Fue.-Air, Turust-Air RATIO 


In contrast to other propulsion devices the com- 
pressor-turbine propulsion relies very little on the 
momentum of the stoichiometric (saturated) fuel-air 
mixture but preponderantly on the much greater total 
mass of air flowing through the duct. 

Assuming, to show this, a work output of 1.5 X 107 
ft.lbs. for each pound of gasoline fuel and a minimum 
requirement of 16 lbs. per sec. air for 1 lb. per sec. 
gasoline, one obtains for the power input Qm’g(m’g = 
p-V.A.g, total weight of air traversing the duct per 
second). 


Q(m'g) (1.5 x 107/ 16) (m’g) combustion 


where (’g) combustion iS the minimum weight of air per 
second necessary for combustion. 
Hence, 


(m 8) combustion = Q x 0.939 X 10-6 (4) 
(m 2) total 


The numerical significance of this result can be esti- 
mated by means of the efficiency formula, Eq. (3), as 
follows. 

The exit pressure ratio g, can be assumed of the 
order of magnitude 1, the efficiency 7 of the order of 
magnitude 0.2, the velocity head can be estimated to 


2 
h = — = — = 0.875 X 10‘ ft. 


the exit velocity ratio, to y, ~ 2.6. This furnishes 
(see Eq. (3)) 


Q = 2(h/n)(y, — 1) = 14 X 10‘ ft.lbs. per Ib. 
Eq. (4) above gives then 


(m’g) combustion /™’Zrota = 14 X 104 X 0.939 X 10-* = 
0.131 


‘and 


totat 0.131/16 = 0.0082 


At any rate, this high dilution (except in the flame 
holder) makes the change of the equivalent molecular 
weight, and with it of the gas constant R, negligible. 


However, the ratio y of the specific heats may change 
noticeably with increasing temperature. 


THE GENERAL FLOW EQUATIONS 


The Bernoulli equation of steady flow for a stream. 
line filament may be written in the form 


d(V?/2g) + vdp = fuds = (f, — f,)vds (5) 


where 
v = 1/pg the specific (per unit of weight) volume 
f = the force intensity (per unit volume) 
fo = the accelerating part of f 
f, = the retarding part of f 
ds = an element of length of the filament of flow 


From the Bernoulli equation, Zeuner has derived a 
second equation by the introduction of the thermo- 
dynamic relations, taking into account the heat added, 
the heat created by the loss of head, and the physical 
properties of gases—that is, the specific heats and the 
equation of state (gas constant, molecular weight, 
etc.). 

This Zeuner equation, as shown below, also contains 
and establishes the function of the concepts of entropy 
and enthalpy (heat content). 

The Zeuner equation may now be developed as 
follows: The first term on the right side of Eq. (5) 
can be expressed by 


favds = dW (5a) 


where dW is the mechanical work per unit of time and 
per unit of weight on the length ds of travel along 
the streamline. 

For the second term f, vds it may be assumed that 
it is the negative work of friction forces (loss of head) 
and that it is entirely transformed into heat dQ,, so 
that 


f,vds = dQ, (5b) 


where dQ, is the heat energy in mechanical units 
(ft.lbs. per Ib.) created on the length ds of travel along 
the streamline. 

If the external heat-input dQ, is added to the heat 
dQ, due to the loss of head, the total heat input will 
increase the internal energy dU and may do work 
pdv so that 


dQ, + dQ, = dU + pdv (6) 


It may now be assumed that the fluid is a gas in which 
the energy-increase dU depends only on the temperature 
increase dT but not on the volume, since this is exactly 
the case in perfect gases and sufficiently exact in 
the actual gases considered here. It follows then 
that 


dU = dT" =d(pr) =d(p)—* 
= = — (6) 
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where 
A = 1/778.26 (B.t.u./ft.lb.) = 1/426.9 (kg.cal.+ 
kg.m.) the reciprocal of the mechanical 
quivalent of heat 
R = R,2d(gp/M,), the gasconstant of a gas mixture 
R, = 1,544 ft.lbs. per °F. = 848 kg.m. per °C., 


the general gas constant 
g = the percentage portion of one gas in a 
mixture of gases 
the molecular weight of one gas in a mixture 
of gases 
Cy/Co = Y, the ratio of specific heats 


The value of y changes somewhat with the tempera- 
ture—namely, from 1.4 to 1.38 between the tempera- 
tures T = 490° and 2,300°F. absolute, while it is 
practically independent of the pressures used in jet 
propulsion. In the later formulas the expression 
y/(y — 1) = wappears. Its value varies between the 
same temperature limits from 3.5 to 3.67. 

With the high dilutions of the fuel in the combustion 
chamber, the resultant reciprocal molecular weight 
2g,/M, due to the combustion differs so little from 
the molecular weight of the air alone that the gas 
constant can be treated as constant in the case of 
turbojet propulsion. 

Inserting the value of dU of Eq. (6’) into Eq. (6) 
one obtains 


d(Q. + Q,) = — + pdv (a) 
Now if Eqs. (6a), (5a), and (5b) are introduced into 
Eq. (5), and it is observed that vdp = d(pv) — pd, 
Eq. (5) assumes the form 
d + d(po) 
2g Cp — 


= + Q) (7) 


It is convenient to bring this equation into a non- 
dimensional form and use an abbreviation for the 
factor d(pv), appearing often, by using the following 
definitions and abbreviations: 


V/Vo = ¢, velocity ratio 
/ Po = q, pressure ratio 
p/‘po = density ratio 
T/T» = 7, temperature ratio 
bo/*/2p0Vo? = H, static over dynamic pressure at the 
flying altitude 


The subscript 0 in the expressions above denotes the 
values in the undisturbed atmosphere at the flying 
altitude. 


(Gp — = 1) = 


As the change of the gas constant after entering the 
duct can be neglected, the equation of state of the gas 
i nondimensional form becomes 


q (8) 
Eq. (7) changes then to 
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d(g? + Hur) = h-d(W + Q) (7a) 


The retarding (negative) work of the friction forces 
has been expressed thermodynamically in the Zeuner 
Eq. (Eqs. (7) and (7a)). It may, in the Bernoulli 
equation, be expressed as the loss of head caused by 
viscous (Poiseuille) and turbulent, momentum 
exchange-shear (Reynolds) stresses. In order to do 
that, one may content oneself with the integration 
of Eq. (5) over the cross sections of a duct section 
since the proposition to express the turbulent and 
viscous friction for an individual streamline filament 
is prohibitively complicated. 

With this simplification of taking the average over 
a cross section, one can write the well-known relation: 


gh 
,uds = , = — ( 
f,vds = dQ g*ds 5c) 


where ¢ is a coefficient depending on the Reynolds 
Number and D, is the hydraulic diameter of the cross 
section. This one-dimensional proposition can be 
justified if annular, not too wide cross sections are 
used. 

For the value of the coefficient ¢ of turbulent friction 
see Jacob and Erk, Stanton and Pannell, and Moody, 
where ¢ is given in function of the Reynolds Number. 
For example, if T = 1,000°F., D, = 1 in., V = 750 
ft. per sec., and the kinematic viscosity » = 0.13 
sq.in. per sec., then the Reynolds Number R, = 
300,000 and the value of ¢ = 0.015. 

The Bernoulli equation, Eq. (5), and the Zeuner 
equation, Eq. (7), with the ratios and abbreviations 
given above, assume the form (/ length of a duct section 
ds’ = ds/l) 


— hW) + Hrd(log q) + = 0 (5d) 
d(g? — hW + Hur — h-Q,) = 0 (7d) 
Observing that ({//4D,)ds’ = h-'dQ, (according to 
Eqs. (5b) and (5c)) and subtracting Eq. (7d) from 


Eq. (5d) one obtains an equation for the entropy S, 
viz.: 
+ QUT = = 
S = Rlog + C 


This functional value of S is valid for any change of 
state—that is, for any condition of p = F(p)[or r = 
F(qg)]. To make Eq. (9) dimensionless, one has to 
introduce a pure uumber that may be called “entropy 
number” = 


> = S/R = log (r°g"') + G (9a) 
The enthalpy (heat content) is defined by 
dI = d (U+ pv) 


and using Eq. (6’), which implies the property of a 
perfect gas 


dI = = RT (10) 
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so that a pure number can be introduced, which may 
be called enthalpy-number J. 


J = 1/povm = I/RT = + (10a) 


Sand I (also = and J) are the entropy and the enthalpy 
functions for any change of state [p = p(p)] of perfect 
gases, 

For any change of state [e.g., 7 = F(q)], the enthalpy 
is expressed by the entropy, so that in the equations 
of flow, Eqs. (5d) and (7d), pressure and temperature 
ratios can be replaced either by the entropy number 
(2) or by the enthalpy number (J) with the result 
that then either entropy number or enthalpy number 
would appear with the velocity ratio ¢ as the only 
variables of state. This representation, however, is 
not chosen for the analysis in this paper.* One reason 
is that in case of leakage or bleeding of heat (dQ 
negative) the flow system is not thermodynamically 
complete and the local entropy value is no measure 
of efficiency. 


CHANGE OF STATE AND CROSS SECTIONS 


The change of state in a flow, in which the influence 
of transverse conduction and radiation of heat is 
either balanced by opposite transverse momentum 
exchange or very small, is expressed by an equation of 
change of state—that is, a relation between two of the 
variables appearing in the equation of state of the fluid. 

The usual expression of a change of state is given by 
p = F(p) sothat the integral { dp/p or the differential 
vdp appearing in the Bernoulli Eq. (5) can be evaluated. 
The change of state in a flow is governed by external 
and internal heat input or loss dQ, and dQ,, respectively, 
and by the variation of cross sections along the duct. 
This can be seen from the fact that, after the change of 
state p = F(p) has been assumed and also after the 
heat and work input have been given, velocity and 
density along a duct follow immediately from Bernoulli 
and Zeuner equations. On the other hand, velocity 
and density along a duct determine the variation of 
cross section by the continuity equation. 

Conversely, the choice of the cross sections, together 
with the heat and work input, determines the corre- 


*In terms of entropy and enthalpy the Bernoulli and Zeuner 
equations appear in the form 


— h-"(W — Q,)] + H[dJ —— dz] = 0 (5e) 
d[g? — h-"(W — Q)] + HdJ = 0 (7e) 
The polytropic change of state [see Eq. (11)] below is in these 
terms given by 
Zi + (A — w) log Jy = 22 + (A — log J2 (11a) 
Here again the conditions of adiabatic, isothermal, constant 
pressure, and constant density are represented by special values 
of X = n/(n — 1). 
n = \ = adiabatic; = 1, = O, isothermal; 
n = 0, = O, constant pressure; n~! = 0, A = 1, constant 
specific volume (constant density). 
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sponding change of state condition (which at supersonic 
velocity can be confined to a shock-wave layer). 

These relations are expressed by the continuity 
condition ApV = const. or 


Agog = = const, (1c) 


If the change of state « = F(qg) is. prescribed, o, g, 
7, and ¢ follow from the solution of the Bernoulli and 
Zeuner Eqs. (5d) and (7d), and A = A(s) (s abscissa 
along duct) follows then from Eq. (1c). 

If, on the other hand, A(s) is prescribed, a relation 
between the velocity ratio gy, the pressure ratio g, and 
the temperature ratio 7 follows from the continuity Eq, 
(1c), and Bernoulli and Zeuner Eqs. (5d) and (7d) then 
determine their functional values in terms of the 
abscissa s along the duct. 

The flow Eqs. (5d) and (7d) are incomplete as long 
as no condition of change of state is introduced. 

For the case where no heat is put in, generated, or 
lost, the condition of adiabatic change of state has to 
be applied—i.e., pv” = const., or 


w= (10') 


For the case with heat added, generated, or lost, a 
condition mathematically convenient and covering a 
wide range of physically appropriate changes is the 
polytropic change of state, pv" = const., or 


= n/(n — 1) (11) 


= const., 


qt —* = const., 


The cases of adiabatic, isothermal, constant pressure, 
and constant density (or constant specific volume) 
change of state can be envisaged as special or limiting 
cases of the polytropic change. _ 

However, it is simpler to treat them by inserting these 
special change-of-state conditions directly into the 
equations of flow, Eqs. (5d) and (7d), as will be seen 
in the case of the “‘perfect’’ jet engine treated as the 
simplest example. 


THE ‘PERFECT’ JET ENGINE 


The analysis of jet engines may be first applied to 
an engine idealized by the assumption that it is per- 
fectly insulated, that the flow is frictionless, and that 
the power of the compressor is equal and opposite to 
the power of the turbine and can be represented by the 
work per unit of time of a force field (f,) acting on the 
streamlines in the duct section of the compressor and 
respectively of the turbine. 

Moreover, here, as well as in the analysis of the 
actual jet engine, the gas constant R and the value of 
pw = y/(y — 1) = 3.5 shall be assumed not to vary, 
since in fact for a temperature range from 500 to 
2,000°F. the value of R changes less than 0.3 per cet 
and the value of u changes from 3.5 to 3.67, as indicated 
in the derivation of Eq. (7). 

Following from the ‘‘ideal’’ conditions, the flow 
through intake and compressor and through turbine 
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and nozzle will be adiabatic, while in the combustion 
chamber a certain freedom is left by the assumption of 
a polytropic change of state that for steady combustion 
will be near a constant pressure condition—i.e., 
n <1, <1. The value of \, together with the 
amount of heat input, determines the cross-section 
change along the combustion chamber, and its effect 
on the propulsion efficiency will also be discussed 

The analysis will give velocities, pressures, tem- 
peratures, and densities everywhere in the unit in terms 
of the data characteristics. 

These data are: 

The static to dynamic pressure ratio 


H = po/"/s poVo? 
The velocity head 
h = Vo?/ 2g 
both in the undisturbed atmosphere and for the velocity 
7, relative to the unit. 
The compression ratio of intake and compressor 
dee = 


The temperature ratio at the transition from the 
combustion (heat) chamber to the turbine 


= 


The (specific) heat input in the heat chamber per 
unit of time divided by the weight of gas passing a 
cross section per unit of time: 


Q, ft.lbs. per Ib. 
The flow Eqs. (5d) and (7d) in the case of adiabatic 


flow (see Eq. (10’)) and for a specific work input W (for 
intake and compressor) are: 


— + — 73) = WwW = 
= (5 adiab.) 


Eqs. (5d) and (7d) become, for polytropic flow (in 
the combustion chamber) with specific heat input Q, 
but no work input, 


ge? — + HX(r, — 1) = 0 
— + Ault, — 71) =h"Q, 


so that 


(5 polytr.) 
(7 polytr.) 


H(u — — 11) = 
— = —Q*a[d/(u — 


Applying these equations for intake and compressor, 
from infinity in front to the entrance to the exit of the 
compressor : 

(a) Intake and Compressor (Adiabatic): 

(12a) 


dec = Tee 


= 1+ — + W*, (13a) 


where 
= W*, = (14a) 


(b) Combustion (Heat) Chamber: 


Jin = Vee =q*, Q,h- = 

Q* Q*,, 
with the abbreviations c, = A/(u — A) and cg = 

u/(u — 2), 

Pen? = = Pec? — Q*, + 
W*, — Q*nc, (16b) 
dea = Qt = Guta Ten” = (17b) 


From Eq. (14b) it can be seen that the requirement of 
the temperature ratio 7* in front of the turbine, the 
compression ratio g* at the entrance of the heat chamber, 
and the amount of heat input Q*, into the heat chamber 
determines the polytropic exponent m = d/(A — 1) 
and with it, by means of the continuity equation, also 
the cross section change along the heat chamber. 
Eq. (14b) also confirms the experience that a higher 
compression ratio g* and a higher specific heat input 
Q*, creates a higher turbine temperature 7* and so, 
because of the limit of heat resistance of the turbine 
material, enforces a limit of compression ratio and heat 
input—that is, a limit of the performance of the jet 
engine. 


(c) Turbine Plus Nozzle (Adiabatic): 


Jen = er Ten = Ter W*,= W*. W*, Pen = ¥e 


Qete* = QutTu” = 
= gen? — W* — — 7) = (18c) 
and by Eq. (15b) 
1+ Hull 1) (196) 


The stepwise analysis of the entire unit applied 
above, shows that the exit velocity ratio y, by Eq. 
(19c) and the product of exit pressure ratio g,, and 
exit temperature ratio 7, (i.e., g,.7.-") by Eq. (18c) 
is determined by prescribing the heat input Q*,, the 
compression ratio g*, and the turbine temperature 
ratio r*. However, the choice of the separate values 
of g, and 7, remains at first arbitrary. It is sometimes 
assumed that the jet acquires at the exit the pressure 
of the undisturbed outside atmosphere. But this is 
not necessarily so, since a pressure excess at the exit 
will outside set itself in equilibrium with the ambient 
atmosphere by momentum transfer (friction) and by 
converging for subsonic and by diverging for super- 
sonic jet velocity. 

From Eq. (19c) and the thrust Eqs. (2a) and (2b), 
it can be seen that the exit velocity ratio y, increases 
with the specific heat input Q*, and decreases with the 
exit temperature ratio 7,, while the thrust also increases 
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with ¢, and decreases with 7,. Moreover, Eq. (18c) 
shows that 7, decreases with an increase of g*. 

In agreement with engineering practice, this proves 
quantitatively that, with fixed turbine temperature 
ratio r* and with fixed specific heat input Q*,, a higher 
compression ratio produces a higher thrust. 


THE CrosS-SECTION DISTRIBUTION OF THE ‘“‘PERFECT”’ 
Jet ENGINE 


After the velocity ratios g and the pressure ratios q 
are determined, the cross-section areas A or their 
ratios A/A, with the exit area are found by means of 
the continuity Eq. (la). As said before, only the 
axial components of ¢ must be used if there is rotation 
of the flow, but in the case treated here, it is assumed 
that in each duct section the rotation is eliminated by 
countervanes so that the tangential components of 
yg are zero at each end of a section. 

Eq. (1a) can, for adiabatic flow in a section, be written 
in the form 


A, = (1c) 
and for polytropic flow 
A; = (1d) 
which for constant pressure (x = 0) becomes 
A, = (le) 


and if in this case mechanical work input or friction is 
excluded (see Eq. (5d)) 


Ge = Yi A, i A; (1f) 


In this way the results of Eqs. (12a) to (19c) of the 
performance of the unit can be used to determine along 
the unit those area ratios without which the per- 
formance values of pressures, densities, temperatures, 
and velocities would not be accomplished. 

One obtains from Eqs. (lc) and (1d) the successive 
area values as follows: 

The flow from infinity to the compressor exit being 
assumed adiabatic, one has: 


[1 + Hu(1—7.)+ 
[1 + Hu(i — + 


where A,, g,, ¢, are the nozzle exit values. This equation 
determines either A,/A,, from a required specific work 
nput W* of the compressor or, vice versa, the specific 
compressor work input from this area ratio and, of 
course, also from the required compression ratio q*. 

Analogously, the area ratio A,,/A,, of entrance and 
exit section of the combustion (heat) chamber is given 
by 


A? (20a) 


1 + Hu(l — + 
(20b) 


1947 


Finally, the area ratio A,,/A, of entrance to turbine 
and exit of nozzle becomes 


[ 1 + Hu(1 — 7.) + Q* 
1+ Hp(1 — + W* — 


Eqs. (20a), (20b), and (20c) show that all areas (A) 
are determined by one of them, preferably by the exit 
area A, and by the value Wp,gV.A, = P, of the com- 
pressor power, prescribed so that it gives feasible 
dimensions and weight of compressor and turbine. 
It seems appropriate to measure W* in terms of Q* 
(say W* = kQ*, k of the order of magnitude 1). 


] (20c) 


RELATIONS BETWEEN THRUST, POWER, SPEED, 
ALTITUDE, AND Exit AREA OF THE “PERFECT” JET 
ENGINE 


The thrust of a jet engine is not uniquely determined 
by the data of performance. Eqs. (18c) and (19c), 
in fact, show that of the three exit variables of tem- 
perature ratio 7,, pressure ratio g,, and velocity ratio 
¢e One remains arbitrary—that is, only expressible in 
terms of the other two exit variables. It is possible, 
however, to complete the conditions by choosing in the 
continuity Eq. (20c) an appropriate area ratio A,,/A, = 
A,,/A, between the exit areas of heat chamber and 
nozzle. 

The choice of this area ratio or what amounts to the 
same of exit pressure g, or temperature ratio r, must 
obviously, all other data unchanged, depend on their 
effect on the magnitude of thrust. 

Different ways to estimate the effectiveness of jet 
propulsion are possible. One method uses the effi- 
ciency »—that is, the ratio of useful to total work of 
Eq. (3): 


1) + H(1 Je") Tee") (3) 


expresses by means of Eqs. (18c) and (19c) two of the 
variables ¢, g., T, by one of them, for instance, ¢, and 
ge by tr, and by means of dn/dr, = 0 determines the 
value of 7, (and ¢, and g,), which gives the maximum 
of 7. 

This method leads to a mathematically involved 
relation between 7, and Q*, which, as in the case of the 
screw propeller, may give the trivial result that the 
highest efficiency appears at lowest power input. 

A second method taken from the theory of the screw 
propeller requires the maximum of thrust Th. if, with 
the other data of speed and altitude, the heat power is 
assumed as prescribed. This method may be repre- 
sented as follows: 

The thrust expression taken from Eq. (2b) is: 


Th. = poVo?A (Ge 1) + (de 1)] (2b), 


The heat-power input P, or better its dimensionless 
parameter II, 
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II = (21a) 


becomes 


P = Q*hgo. VA. 


Tl = Q*qepere™ (21b) 


One has then to determine the maximum of thrust for 
a prescribed power by requiring the value zero for the 
differential coefficient after g, and y, are expressed by 


Te: 


bs [Th. + LM] = 0 (22) 
dr, 


where L is the Lagrange multiplicator to be determined 
by the above Eq. (21b) for II and after g, and ¢, have 
been expressed by 7, by means of Eqs. (18c) and (19c). 

The author has not yet succeeded in deriving from 
this second proposition sufficiently simple conclusions. 

A third proposition requires the thrust coefficient 
Th.* (that is, the thrust per unit of nozzle exit area 
A, and of ambient dynamic pressure '/2poVo”) to 
become a maximum for a prescribed specific heat 
input Q* = Qh-'—that is, for a heat power input per 
unit of weight flow p.gV.A, and per unit of head h = 
Vo?/2g of the ambient atmosphere. Both quantities 
Th.* and Q* are dimensionless parameters. 

This last proposition, which leads to fairly simple 
results, is in detail illustrated by an example. 

A thrust coefficient Th.* = Th./po Vo?A, was given in 
Eq. (2a) in the form 


Th.* = — 1) + (qe —1) (2a) 


The results of the performance Egs. (12a) to (19c) 
serve to express the exit temperature ratio 7, and the 
exit jet-velocity ratio g, by the exit pressure ratio q,. 
In these expressions the coefficients depend on the 
prescribed values of Vo, po, po, Zo; on the prescribed 
compression ratio g*; on the prescribed turbine tem- 
perature ratio r*, and on the specific heat-input Q*. 

The relations between g,, ¢,, and 7, may be repeated 
in the form 


ge = (18¢) 
where 
= k — Hur, (19¢) 
‘where (see Eq. (19)) 
k=Q*+1+ Hp (19’c) 


Eq. (2a) can then be expressed in a form that con- 
tains as the only variable r,—namely, 
Th.* = — Hur, — (k — Hur.)'*] + 
(2b) 


To obtain an extremum value for Th.* one must 
make vanish the differential coefficient with respect to 
7, Which leads to: 
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dTh.* 
dr, 
Canceling the factor r,"~—', which is certainly different 
from zero, and multiplying by (k — Hur.)'", which, 
according to Eq. (19’c), is also not equal to zero but is 
greater than 1, one obtains the equation 


(k — Hur, — 1)[Hu(2u — — 2k(u — = 0 


or, since the first factor (see Eq. (18c)) has the value 
Ge" i> 0, 


= — — Hur, — (k — + 


=x — = 2(u — 1) 
(23) 


From Eq. (19c) the exit velocity ratio ¢, is then given 
by 


and from Eq. (18c) the exit pressure ratio by 
qe = (25) 


The maximum thrust coefficient then becomes 


{1 — — 1) (2c) 


However, since the expressions for the exit values of 
temperature, relative jet velocity, and pressure ratio are 
given by Eqs. (23), (24), and (25), it seems to be more 
convenient to use, as will be done below, the form of 
Eq. (2b) instead of Eq. (2c) for the thrust coefficient. 

The value of the polytropic characteristic 4 = 
n/(n — 1) (n polytropic exponent) in the combustion 
chamber is in general assumed as zero, which means 
constant pressure change of state. The question is 
whether this value of is the best one for the thrust 
maximum. Looking at the expression for 7,, ¢,, and 
ge to find an answer to this question, one sees that 
only g, is dependent on the value of \ (see Eq. (25)) 
because of the term] = (g*~/“r*)!-°/™, One also 
sees that \ > 0, m < 0 gives a smaller value of /~ and 
therefore a smaller value of Th.*, so that Eq. (2b) 
(and also Eq. (2c)) give smaller thrust values if \ 
is positive—that is, tending toward p = const. 

On the other hand, a doubt about the allowable 
value of excess of exit pressure ratio g, does arise, since 
it must not be higher than the pressure excess realizable 
by an equilibrium between free jet and ambient air 
arising from friction (i.c., momentum transfer) and 
from convergence of a subsonic or divergence of a super- 
sonic free jet. 

The case \ < 0, > 0, which deviates from constant 
pressure toward isothermal change of state, seems to 
give increase of thrust but has not yet been evaluated. 
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The thrust itself is obtained from Eq. (2b), which so that 


may be repeated here 
(Th.)maz- = (Th.*) maz,poVo7A. (2a) 


Numerical values of thrust coefficient and thrust are 
now ready for calculation from the data. Using for 
instance the data given below, one has to proceed as 
follows: 

Compression ratio and turbine temperature ratio 
may be required to be g* = 4and7* = 4. The change 
of state in the heat chamber may be required to be of 
constant pressure—that is, n = 0,’ = 0. Witha 
velocity Vo = 750 ft. per sec. and an altitude of 22,000 
ft., 


Po Vo? 
= 2.82, h=— = 0.875 X 10! 
1/299 Vo? 2g 


Hy = 2.82 X 3.5 = 9.87, (I),-0 = 4074 = 2.7 
For the specific heat input (see Eq. (15b)) 
Q* = H(p — d)(r* — = 9.87(4 — = 24.8 


It follows then 
k = Q*+ 1+ Mp = 35.67 


k (*—— ) = 35.67(—— ) = 5.94: 


¥=4 ) 71.34/ 0A ) 
7 = = = = 3.01 
H\y(y + 1) 2.82\1.4 X 2.4 


de = — = = 1.558 
14 X 2.82: X. 2.7 


The thrust coefficient with these values (see Eq. (2b)) 
becomes 


Ge” 


1.558 .945 — 2.435 


(1.41 X 0.558) = 1.811 + 0.786 = 2.597 


Assuming a nozzle exit-area A, = 1.5 sq.ft., the thrust 
will have the value 


Th. = 2.597p)Vo?A, = 2.597 X 0.0012 XK 750? XK 1.5 = 
2,620 Ibs. 


If now the thrust value corresponding to the exit 
pressure ratio g, = 1 (p, = po) shall be determined, in 
order to compare it for the same data with the maximal 
value calculated above, one has to proceed as follows: 

The equation for the thrust coefficient is simplified to 


From Eq. (18c) can be seen 


From Eq. (19c) 
= Q* + 1+ — = 9.01 


(Th.*)gem1 = (1/2.7)6.01 = 2.225 


(Th.)gem1 = 2.225 X 0.0012 XK 750? X 1.5 = 2,245 
Ibs. The gain in thrust is for the specific heat input 


of 
Q* = 24.8 


(Th.*) max./(Th.*)ae=1 (Th.) maz./(Th.)ae=1 = 1.168 


From Eq. (2b) and (3) the efficiency formula can be 
stated in the form 


With the values computed above one obtains: for 
ge > 1, n = 16.65 per cent; for g = 1, 7 = 16.15 
per cent. 

The comparison of the efficiencies n for the case of 
equal specific heat input Q* but for different exit pres- 
sure ratios g, shows the peculiar result of smaller 
efficiency at smaller exit pressure but at still smaller 
thrust. 

The reason is that in the efficiency expression (see 
Eq. (21b)) the value Q*q,¢,7,~—' of the total heat-input 
power appears in the denominator. This may cause, 
with increasing exit pressure ratio q,, the factor g,¢,7,~1, 
denoting the mass flow, to increase more than the 
thrust-power value in the numerator. 

This observation, together with the consideration of 
the feasible values of exit pressure augmentation by 
friction, transverse momentum transfer, and transverse 
centrifugal force of the converging streamlines of the 
jet emerging into the ambient atmosphere, indicates 
a certain advantage for maximum method 2 because it 
introduces the total heat power input (see Eq. (22)). 


ANALYSIS OF LOSSES IN THE JET ENGINE 


The Nonadiabatic Bernoulli-Zeuner Equations of Flow 


The perfect (or ideal) jet engine was analyzed by 
means of Eqs. (5d) and (7d) with the simplification of 
omitting the turbulent loss of head ({y?/4D,) per unit 
of streamline length, the actual excess of power of 
turbine over the power of compressor, and the losses 
of heat due to imperfect insulation. 

These facts, making the flow in all sections of the 
duct unit nonadiabatic, necessitate an assumption for 
the changes of state in the different sections of the unit, 
which afterwards must be justified by the energy 
balance and the choice of the cross sections along the 
duct. 

Following this reasoning, a polytropic change of 
state 

gr~* = const. (11) 


with possibly different in each section of the unit, 
may be introduced. 
Eq. (10) is equivalent to 


= (1 1b) 
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where the subscripts 2 and 1 denote two cross sections 
of one section of the unit. 

It is also well known that all of these relations have 
the same structure as the adiabatic relations inasmuch 
as the difference between them consists only in changing 
yinto mor wp = y/(y — 1) intod = n/(n — 1). 

Using Eq. (10b) to eliminate g from Eq. (5d), one 
obtains the equations 


d dw 
as! (g? + + 2¢ as’ (5e) 


d d 
— (¢? + Hur) = — (W* + Q*.) (7e) 
ds ds 

where, for abbreviation (see Eq. (5c)), 2 = §1/4D,. 


The Integration of the Bernoulli-Zeuner Equation 

The mechanical (positive or negative) specific work 
input W*, heat input Q*, and coefficient 2(s’) may be 
assumed to be prescribed. Then Eqs. (5e) and (7e) 
appear as two equations with the two unknown 
variables +r and yg. The temperature ratio function 
r(s’), eliminated from the two flow equations above, 
leaves, for the square of the velocity ratio g, the 
differential equation 


de" 

ds’ 
With the abbreviations 

Ce — = A/(u — 
the integral will be 


+ = (Wt -aQ%) 


(d/ds’) (W* — ors | (5g) 


where the subscript 7 denotes the entrance station of 
the duct section in consideration. If one is content 
with the determination of the state at the entrance 
and exit of each duct section, the integral functions 
become constants. It is convenient to introduce the 
following abbreviations: 


b = ((W* — — 
(27) 


Eq. (5g) specialized for the values at the entrance (7) 
and the exit (e) station of any duct section then assumes 
the form 


= + (W* — GQ*)ob (5h) 


Eqs. (Se) and (7e), after subtraction and in consequence 
of Eq. (5h), furnish 


Hur, = Hur, + 92(1 — a) + [(W*(1 — 6) + 
Q*(1 + (7h) 
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In the following development W* and Q* may denote, 
respectively, the total specific (positive or negative) 
work and heat input of a section of the unit, whereby 
the different sections may be indicated by the sub- 
scripts i for intake, c for compressor, h for heat chamber, 
t for turbine, and for nozzle. 7 

The analysis can now be performed step by step 
analogously to the analysis of the perfect jet engine. 
Only the number of steps will be higher, since for 
nonadiabatic flow neither intake and compressor nor 
turbine and nozzle can be treated as one step. To 
be quite exact the heat chamber should also be treated 
in the two sections of flame holder and heat flow, but 
this must be left for a special investigation. However, 
if only the states at the exit and entrance are to be 
considered, the treatment of the heat chamber in one 
step, with an average change of state, may be justified 
as a first approximation. 

To the two last flow equations, Eqs. (5h) and (7h), 
the change of state condition between the entrance and 
the exit stations of each duct section must be added as 
the third equation—namely, 


qi(te/71)” (11b) 


The continuity Eq. (lc), combined with the pre- 
ceding equation, provides a fourth relation between 
the entrance and exit cross-section ratio and the 
pressure and temperature ratios at these same cross 
sections—namely, 


= (1d) 


where the polytropic exponent m for adiabatic change of 
state becomes y, for constant pressure » = 0, for con- 
stant temperature n = 1. 


SUCCESSIVE ANALYSIS OF A ‘‘NONPERFECT”’ JET ENGINE 


The four equations, Eqs. (5h), (7h), (11b), and (1d), 
are consistently arranged in such a manner that the 
significant variables at the exit (e) of a duct section 
follow immediately from their values at the entrance 
(i). The exit values then serve as the entrance values 
for the adjoining duct section. 

Thus, the analysis is as follows: 

(1) The Intake.—First subscript 7. Since the flow 
arrives from the undisturbed atmosphere (i.e., g, = 1, 
t~.=1,¢,, = 1) and since no mechanical work or external 
heat W and Q is put in, Eq. (5h) furnishes (first sub- 
script c for compressor) 


Pie? == Pie? = 
a, = S0' sds’ (Eq. (26)) (8) 
Hyuty, = Hyt., = Hu + 1 — a; (29) 
die = = Tie = Tea +1-— a;)/Hp)™ 
(30) 
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There is a certain interest in providing, at the 
entrance of the compressor, a high pressure ratio 
die = In order to investigate this possibility, 
a numerical example is given. 

Assume again Hu = 9.87, Vo = 750 ft. per sec., 
po = 0.0012, ¢ = 0.015, 1/D, = 10. Then one obtains 


a= 


Some preliminary estimates show that appropriate 
values of the polytropic exponent are 


3) = 1.5 
nm, = 5, \y = 1.25 


The table below shows the consequences of these 
assumptions. 

m Cus ai Pie tie Aie/Ais Yin 

3 1.5 1.75 0.013 0.114 1.10 8.4 1.153 

6 1.25 1.556 0.02 0.1414 1.0995 6.9 1.126 

The result of the last column of the table above 
following from the equation of continuity, can be 
interpreted either envisaging A,, as the cross section 
in front of the intake or as the cross section of the 
incoming free column of air far in front, where the 
disturbance of the atmosphere by the jet engine is 
negligible. 

The first interpretation, if required with appreciable 
pressure increases, would lead to flow separation and 
backflow caused by the corresponding diameter increase 
at feasible intake length. 

At the front station of the intake in the latter inter- 
pretation, there is already a partial stagnation so that 
the increase of areas along the intake proper is smaller 
than appearing in this last column of the table. With 
this interpretation the front area at infinity may be 
called Ao. 

It is evident that by prescribing the ratio of diameter 
increase to intake length one can calculate velocity 
ratio g,; and pressure ratio q;; at the front section 
of the intake. 

However, for a given exit area of the intake (i.e., 
entrance section of the compressor) the resultant mass 
flow of air, which is important for the attainment of 
thrust, becomes smaller if the intake and the counter- 
pressure of the compressor is designed for higher stag- 
nation (ram) pressure. 

(2) The Compressor.—The heat that the compressor 
produces by compression and by friction leaves again 
partly through the inner and outer walls. It is plausible 
to assume for the loss of heat a percentage of the specific 
work (per pound weight and per unit of length of a 
streamline) of the compressor, which, with the sign 
definition of this paper, is positive. This loss, there- 
fore, is expressed by 


Q*, —6,W*, 
so that 
(W*, = W*d.(1 + 


and 
W*.(1 — + Q*-(1 + 4.6.) = W.[1 — 6(1 + 


If now for convenient writing one introduces the 
abbreviations 


(1 —a,)/Hue =a, + = 8, 
and 
(1 — 6. — &)/Hu = 
Eqs. (5h) and (7h) become 
Gee? = + W*.B, (32) 
Tee = Tei + + W* (33) 
The change of state Eq. (0b) gives 
Tee Gee = Tei ™ Get (34) 
and the continuity equation 
dee!” GeeAce = hes (35) 


Using now the values of q.;, ¢:% 7-1, and A,, found 
for the intake, one obtains by means of Eqs. (32) to 
(35) the exit values of velocity, temperature, pressure, 
and either area ratio A,,/A,; by the polytropic ex- 
ponent , or the polytropic exponent m, by a chosen 
area ratio A,,/A,, and the specific work-input 
W*, of a value limited by permissible weight and 
length. 

If, however, as was done for the ‘‘perfect”’ jet engine, 
the compression ratio g* is prescribed, then all four 
variables ¢,,, Te, Ace, and m, are expressed in terms 
of 

(3) The Heat Chamber.—In the heat chamber no 
work is performed (W*, = 0), but the heat input 
Q,h-! = Q*, may be diminished by insulation losses 
proportional to it of value —6,Q*,, so that 


Cyn bn(Q* neftect. = — = 

where Q*, is the original heat input and 

Bn = — 
Furthermore 
O* neftect.(1 + = O*n(1 — 6,)(1 + Cyngn) = 
where 

ya = (1 — 6, + By)/Hu 

and 

a, = (1 — a,)/Hp 


The flow equations then become 


Pre? = Q* (36) 
the = Tri + + O* arn (37) 
The *Qne = Tri OF Gy, = (38) 
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(32) 
(33) 


(34) 


These equations again determine velocity, temperature, 
and pressure ratios and either polytropic exponent or 
area ratio in terms of the corresponding determined 
values of the preceding section. If, moreover, as for 
the “‘perfect’’ jet engine, the temperature ratio r,, = 7* 
at the exit of the heat chamber is prescribed, all four 
variables have fixed values only depending on one of 
the areas. Eq. (37) analogously to Eq. (15b) de- 
termines again the heat input in terms of compression 
ratio and turbine temperature. 

(4) The Turbine.—In the turbine no heat is added, 
but heat of amount proportional to work produced 
(W, negative) expressed by Q*, = W*,6, may leave. 
If, then, for abbreviation one introduces 


and 
(1 + 6, — = 
a, = (1 — a,)/Hp 


the flow equations assume the form: 


Pre” = + W*,B, (40) 

Tre = Tre + + (41) 
= QneTre™* (42) 
™ = Qne/™ PncA ne (43) 


One has again five new quantities [i.e., g:., Tre, Qres 
A,, or \; = n,/(n, — 1), and W*,] but only four equa- 
tions. However, there is one more relation between 
the specific work input W*,of the compressor and the 
work output W*, of the turbine—namely, 


W*, W*, € (44) 


where «€ > 1 expresses the additional work of the 
turbine to be performed for accessory motors or gener- 
ators, for the covering of friction losses in bearings 
and gears, and, eventually, for a screw propeller. 

With this relation the turbine variables are all 
expressed in terms of the corresponding values of the 
preceding (heat chamber) section and of the specific 
compressor work. 

(5) The Nozzle.—In this section neither work nor 
heat is put in. Only a heat loss may appear, which 
here may be assumed proportional to the exit tem- 
perature-ratio 7,, of the turbine expressed by 


= —Trdp 
The flow equations are now 
= Ge? = (45) 
Tre =Te = Tre + — Tren (46) 


Where y, = 6,(1 + and a, = (1 Qn) / 


Te "Ge Tre (47) 
de = gre (48) 


ANALYSIS OF THERMAL TURBOJET PROPULSION 


One has again here the five unknown quantities ¢,, 
Te, Jer Xn = Mn/(M, — 1), and A, and only four equations 
to determine them. 

However, as in the analysis of the “‘perfect’’ jet 
engine, an additional relation between g, and 7, must 
be prescribed. This must be done either by a maximum 
condition for the thrust at prescribed total heat input 
or by the condition of the exit pressure equal to the 
undisturbed atmospheric pressure—that is, by the 
condition g, = 1. 

If this is done, then all state variables (velocity, 
temperature, pressure, and area ratios and polytropic 
exponents) are expressed in terms of either the intake 
entrance A, or the nozzle exit area A,; the specific 
compressor work |W*,; the specific heat input Q*, in 
the heat chamber;. the compression ratio g* = q..; 
the turbine temperature ratio 7* = 1,;; also the data 
HT, h, R, and y and the section-wise loss coefficients a, 
b, 6, and e; the gas constant, R; the adiabatic ex- 
ponent y = c,/c,; and the losses due to leakage and 
internal drag expressed by the coefficients 8, 4, 
and Yicnin- 
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Stresses in Cylindrical Semimonocoque 


Open Beams 


HENRY L. LANGHAAR* anp CLARENCE R. SMITHT 


Consolidated Vultee Aircraft Corporation 


SUMMARY 


The theory of open, symmetrical, semimonocoque beams of 
circular-cylindrical shape is developed by direct minimizing of 
the strain energy in a length of open beam whose ends are acted 
upon by specified resultant forces and moments. Simple explicit 
formulas are obtained for the stresses. The theory exhibits 
exponential stress gradients near the longerons—an effect that 
has not heretofore been incorporated in the analyses of open 
beams. This phenomenon may cause high axial stresses at the 
corners of a cutout in a cylindrical beam and high shearing 
stresses near the longerons. Experimental data are presented 
which tend to confirm these conclusions. 


INTRODUCTION 


Sipe THEORY of open semimonocoque beams is of 
interest to aircraft designers because it affords a 
basis for the stress analysis of structures such as the 
portion of a fuselage in the interval of a cargo door. 
The experimental data presented in conjunction with 
the following theory of open beams indicate that this 
theory is satisfactory for these applications, even 
though the cutout is short. The transition from the 
open-beam distribution to the closed-beam distribution 
has not been theoretically investigated in this study. 
However, the test data indicate that the transition 
length is fairly short, even if the beltframes are ex- 
tremely flexible. This suggests that hoop stresses are 
active agents for effecting redistribution in the region 
beyond the cutout. 


NOTATIONS 


Consider a symmetrical open beam of circular cross section. 
as shown in Fig. 1. Let the x, y, z axes form a right-handed 
system (positive z-axis directed toward the reader). (Fig. 1.) 
Let the plane z = 0 lie at the mid-section of the open beam. Let 
the beam be subjected to a vertical shear V through the geometric 
axis of the cylinder and a torque M (Fig. 1). Let the resultant 
x- and z-components of the force and the resultant y-component of 
the bending moment be zero. Let the resultant forces, bending 
moment, and torque at any station 2 be the resultants of the 
actions of the portion of the beam in the region z > 2; upon the 
section at 2. The x-component of the bending moment at 
station z, then has the value + 

The following nomenclature is employed: 


= skin thickness (assumed to be constant) 
equivalent thickness of the plate-stringer combination 
(assumed to be constant). The equivalent thickness 
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may be defined by the formula ¢, = ¢ + (A,/b) in 
which A,/b is the ratio of the cross-section area of a 
stringer to the stringer spacing 

cross-section area of either bar (longeron) reinforcing 


a = 
the free edge (Fig. 1) 
2L = length of the open beam 
r = radius of the cylinder 
28 = central angle included by the cross section of the beam 
(Fig. 1) 
6 = polar coordinate (Fig. 1) 
N_ = axial load per unit length of the plate-stringer cross- 
sectional periphery (tension positive) 
q = shear load per unit length of the plate-stringer cross- 
sectional periphery (positive counterclockwise, 
Fig. 1) 
o = N/t, = axial stress in the skin and stringers 
t = g/t = shearing stress in the skin 
Mzo = bending moment about the x-axis at the station ¥ 
z = 0 (positive when it tends to rotate the section 
in the sense of a right-hand screw advancing along 
the x-axis) 
= shear load in the direction of the positive y-axis 
(Fig. 1) 
M = torque with reference to the geometric axis of the ; 
beam (positive counterclockwise, Fig. 1) 
No = value of N at the station z = 0 
T = tension in the lower longeron 
v = Poisson’s ratio 
E = Young’s modulus 
G = 0.5E/(1 + v) = shear modulus 
U- = strain energy in the open beam 
d = a/rte = a dimensionless parameter 
k? = 6r%, (1 + v)/L% = a dimensionless parameter 
Izz = 1%t,(8 — sin B cos B + 2d sin? 8) = moment of inertia 
of the beam cross section about the x-axis 
4 
| 
>-X 
6 
V 
Fic. 1. Cross section of open beam showing nomenclature. 
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STATICAL CONDITIONS 


It is assumed, as is usual in open-beam theory (see 
references 1 and 2) that hoop stresses in the skin are 
negligible. The differential conditions of equilibrium 
for a rectangular element then show that the shear q 
is independent of the axial coordinate z and that the 
axial force JN is a linear function of 2; specifically 


= (—2q'/r) + No (1) 


in which the prime denotes a derivative with respect 
to the angular coordinate @. 

The differential condition of equilibrium for the 
lower longeron (Fig. 1) is dT /dz = q, in which T is 
the tension in the longeron and gq is the shear in the 
adjacent skin. Also, continuity requires the relation- 
ship 


o = T/a = N/t, (2) 
in which JN is the axial load intensity in the plate- 


stringer combination adjacent to the longeron. Hence, 
by Eq. (1), 
dT/dz = —aq'/rt, = — dq’ 
in which 
h = a/rt, (3) 


Thus, the boundary condition at the lower longeron is 
q+ rq’ = 0 (4) 

The corresponding condition for the upper longeron is 
q— dq’ =0 (4’) 


In addition to the preceding differential conditions 
of statical equilibrium, there are six integral relation- 
ships of statics. The first three of these relationships 
pertain to the shears, and, since shear loads in the 
longerons may be assumed to be negligible, the equa- 


tions do not involve the longeron forces. Hence, the 
equations are: 
—V/r 
S',q sin 0d0=0 (5) 


qd0 = M/r? 


In setting up the equilibrium conditions for the 
axial forces, the tensions in the upper and lower 
longerons at the mid-section are respectively designated 
by —7) and Ty. This is in conformity with an anti- 
symmetry of the plate-stringer stresses which is subse- 
quently confirmed. It follows: 


Nod 0 = 
2(T)/r) sin B ++ No sind dé = 
No cos 6d 6 =0 


—My/r’ (6) 


STRAIN ENERGY IN AN OPEN BEAM 


The strain energy in the beam is represented by the 
equation 
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With Eq. (1), this may be expressed 


in which an irrelevant constant factor is discarded. 
This may be integrated with respect to z, and there 
results 


U= 


2(1+7) | x 
d 0 + K U 


2(1 + v)(t.q?/t)] 
(7) 


+ Not + 


in which the factor 2Z is discarded. 


MINIMIZING OF STRAIN ENERGY 


If the strain-energy expression, Eq. (7) is minimized 
with respect to No, under the restrictions imposed by 
the auxiliary equations, Eq. (6), the elementary flexure 
formula results;* i.e., the axial stresses at the mid- 
section of the open beam (i.e., the section z = 0) are 
expressed by the ratio of the bending moment to the 
section modulus. This result may be expressed ana- 
lytically as follows 


No = —(Myort,/I,:) sin (8) 


It is convenient to introduce a constant y, defined by 
the equation 


¥ = M,/[r?(6 — sin Bcos B + 2dsin? B)] (9) 
Then 
No = —ysin 0 (10) 


To minimize the strain energy with respect to q, an 
undetermined linear combination of the integrals of 
Eq. (5) is added to the strain-energy expression, Eq. (7) 
This is a common adaptation of Lagrange’s method* 
of undetermined multipliers for the treatment of 
variational problems of this type. The modified 
strain-energy expression is 


= [(L2q’2/3r2) + 2(1 + v)(t,92/t) + org cos 6 + 


in which 4, ¢2, 3; are constants. 

The minimizing of U is a variable end-point problem 
in the calculus of variations. According to Eqs. (4) 
and (4’), variations of q are restricted by the end 
conditions. 


+ = — = O 


However, it may be shown that these conditions are 
necessary for a relative minimum of U when q is freely 
varied. This result might have been anticipated; 
it merely verifies that Eqs. (4) and (4’) could have 


* See Appendix. 
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been derived as natural boundary conditions of the 
variational problem. 

Irrespective of the end conditions, the Euler equation 
for g must be satisfied. This equation is 


2L°q"/3r? — 4(i + v)(t.g/t) — o1 cos 
— = 0 (12) 
Let 
k? = Gr*t.[(1 + v)/L*t] (13) 
Then, with a modification of the ¢’s, Eq. (12) may be 
expressed : 


— = cos 6 + sin + (14) 
The general solution of Eq. (14) is 


q = B, cosh ké + sinh + Az cos 6 + A; + A, sin 
(15) 


in which the A’s and B’s are constants. 

The second of Eqs. (5) is automatically satisfied if q 
is an even function of 6. Also, this condition satisfies 
the boundary condition for the upper longeron when 
Eq. (4) is satisfied. Therefore, a solution will be sought 


with the assumption that g is an even function. Then 
B, = A, = 0, and Eq. (15) reduces to 
q = B, cosh k@ + Az cos 0 + A; (16) 
The first and the third of Eqs. (5) then yield 
2(1 + k?)~'(sin cosh kB + B sinh k8)B, + 
(8 + sin 8 cos 8)Ay + 2A3sin 8B = —V/r > (17) 
B,k— sinh k8 + Az sin B + A38 = M 


The boundary condition Eq. (4) yields 


B,(cosh kB + Ak sinh kB) + A2(cos B — Asin B) + 
A; = 0 (18) 


Now sinh k8 ~ cosh k8 0.5exp It is there- 
fore convenient to introduce a new constant defined by 
A, = 05B, exp k8. Then Eqs. (17) and (18) become 

2(1 + k*)—(sin 8 + k cos B)Ai + (8 + 

sin 6 cos + 2A3sin B = —V/r 

Ao sin B + = M/2r? 

A,(1 + AR) + A2(cos BA sin B) + As = 0 


(19) 


The first term in Eq. (16) is negligible when @ is 
small. Hence, it is permissible to make the approxi- 
mation 


cosh k@ = 0.5 exp 


in which the symbol |@| denotes the absolute value 
of 6. Then, Eq. (16) becomes 


q = 4 A, cos 0+ As (20) 


in which e is the base of natural logarithms. Hence, 
with the aid of Eqs. (1) and (10), the axial load V may 
be expressed by the equation 


N= (sr~'Ay — sin @ (21) 


in which the factor (+1) is opposite to the sign 
of 0. 

Eqs. (19), (20), and (21) are a solution of the problem. 
The A’s are determined by the linear algebraic equa- 
tions, Eq. (19). These constants are to be substituted 
in Eqs. (20) and (21), to furnish the shears and the 


Fic. 2. Photograph of test setup. 
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axial loads. The stresses are determined by the equa- 
tions r = g/tando = N/t,. 
vont 
TEST SPECIMEN AND PROCEDURE me STA7 S000LB 
In order to test the foregoing theory, a 24S-T Alclad 
cylinder 30 in. in diameter and 168 in. long was con- 10000 ra (el Pcie 
structed. The skin was 0.027 in. thick. The cylinder sooo|_|_| fal u 
was reinforced by internal extruded angle-type stringers o STA.6 
of cross-section area 0.087 sq.in. spaced at 15°, and S 
external strap-type 1/, in. by 15/s in. hoops spaced at td Rel 
14 in. 5000 aN 
The cylinder was set up as a cantilever beam. Shear ¢ o S7A. 
loads were applied by means of a calibrated hydraulic 5 s / a 
ram. Fig. 2 shows the test setup for the cantilever = aN 
tests. An adaptation for making torsion tests was . 5 NL 4 of ‘ 
also installed. This was accomplished by radially > SAF 
displacing the loading ram 15 in. and reacting the . 
applied load by a symmetrically located turnbuckle My, 
in series with a spring dynamometer. (Fig.2.) Thus, a. eT ; 
a couple with a 30-in. moment arm was applied to the pe 


cylinder. Fig. 3 is an elevation view of the cylinder 
showing station locations and the point of application 
of the shear load. 


+ 
SERES 


_| 
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Fic. 3. Elevation drawing of cylinder, showing 30° by 14-in. 
opening. Dotted lines indicate subsequent enlargements of 
opening to 60° by 14-in., 90° by 14-in., and 90° by 28-in. 


Stresses in the test cylinder were determined by 120- 
ohm electric strain gages. Strain gage rosettes on 
opposite sides of the skin were series-connected, in 
order to eliminate bending effects. Single gages were 
attached to the stringers to measure axial stresses. 

A test run consisted of two sets of readings of the 
strain gages. These readings served to determine the 
slopes of the load-stress curves. In all cases, the loads 
were controlled so that the proportional limit of the 
material was not exceeded. 

In the first test, the cylinder had no opening. This 
test served as a calibration of the strain gages. After 
the first test, a 30° by 14-in. opening was cut near the 
center of the cylinder. This opening was framed by 
two uncut stringers and two consecutive hoops. In 
the course of the tests, the opening was enlarged and 
various types of edge-reinforcements were installed. 
The following is an outline of the test conditions: 


Test No. 1—No opening. 
Test No. 2—30° by 14-in. opening bordered by uncut 
stringers and hoops. 


-60 60 120 180 


Fic. 4. Stringer stresses for Test No. 1 with cantilever loading, 
Curves are derived from the elementary beam formula, ¢ = My/I, 


Test No. 3—60° by 14-in. opening bordered by uncut 
stringers and hoops. 

Test No. 4—90° by 14-in. opening bordered by uncut 
stringers and hoops. 

Test No. 5—Same as Test No. 4, but with '/, in. by 
in. 248-T  edge-reinforcing strips 
(longerons) extending one bay beyond 
each end of the cutout. 

Test No. 6—Same as Test No. 5, but with 1/, in. by 
15/s in. steel hoops replacing the 24S-T 
hoops at the ends of the opening. 

Test No. 7—Same as Test No. 4, but with 1/, in. by 
15/; in. steel hoops replacing the 24S-T 
hoops at the ends of the opening. 

Test No. 8—90° by 28-in. opening bordered by uncut 
stringers and hoops. The hoop at the 
center of the opening was cut out in the 
region of the opening. The center hoop 
and the two hoops bordering the opening 
were 1/, in. by 15/s in. steel straps. 

Test No. 9—Same as Test No. 8, but with 24S-T 
1/, in. by 15/s in. hoops replacing the 
steel hoops. 

Test No. 10—Same as Test No. 9, but with 1/, in. by 
in. 248-T edge-reinforcing strips 
(longerons) extending one bay beyond 
each end of the opening. 


DISCUSSION OF RESULTS 


Theoretical stress curves, determined by Eqs. (19), 
(20), and (21) are shown in Figs. 4 to 17. Test points 
are plotted on these figures. Test points for the bays 
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Fic. 5. Stringer stresses for Test No. 2 with cantilever loading. 
Curves for Stations 5 and 6 are derived from Eq. (21). 
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Fic.6. Stringer stresses for Test No. 3 with cantilever loading. 
Curves for Stations 5 and 6 are derived from Eq. (21). For 
comparison, a curve derived from the formula ¢ = My/I is also 
shown for Station 6. 


beyond the cutout are also shown, and the theoretical 
stress curves for the closed beam with no opening are 
plotted for these bays in order to indicate the rate of 
transition to the undistorted stress pattern. 

The test points exhibit considerable scatter. This 
may be partly attributable to incomplete elimination 
of local bending stresses. The fact that Fig. 4 shows 
two stringers that gave consistently low stresses at 
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Fic. 7. Stringer stresses for Test Nos. 4 and 7 with cantilever 
loading. Curves for Stations 5 and 6 are derived from Eq. (21). 
Circles and triangles show data from Test Nos. 4 and 7, re- 
spectively. 
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Fic. 8. Stringer stresses for Test Nos. 5 and 6 with cantilever 
loading. Curves for Stations 5 and 6 are derived from Eq. (21). 
Circles and triangles show data from Test Nos. 5 and 6, re- 
spectively. 


four stations indicates that these stringers were bowing 
under the application of the load, and that the bending 
stresses were detracting from the average stresses. 
Despite the test scatter, the data show a marked 
trend to follow the theoretical curves. It is especially 
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= 


- -120 -60 O 
Fic. 9. Stringer stresses for Test Nos. 8 and 9 with cantilever 


loading. Curves for Stations 4 and 6 are derived from Eq. (21). 
Circles and triangles show data for Test Nos. 8 and 9, respectively. 
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Fic. 10. Stringer stresses for Test No. 10 with cantilever 
loading. Curves for Stations 4 and 6 are derived from Eq. (21). 


noteworthy that the tests exhibit extremely high 
stringer stresses at the corners of the opening and high 
shearing stresses near the longitudinal edges of the 
opening. This is in agreement with the exponential 
stress gradients that the theory predicts. 

Fig. 18 shows test points for a run in which the beam 
was subjected to a symmetrical cantilever load that 
produced equal tensions in the two stringers bordering 
the opening. The curves on this figure represent the 
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Fic. 11. Shearing stresses for Test No. 1 with torque loading. 
Curves are derived from the Prandtl torsion formula, r = //2At. 
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Fic. 12. Shearing stresses for Test No. 2 with torque loading. 
Curve for Station 5.5 is derived from Eq. (20). 


theoretical stringer stresses as determined by the 
elementary flexure formula; i.e., the ratio of the bend- 
ing moment to the section modulus. It is apparent 
that the stresses at stations numbers 4 and 7 are not 
fully developed in the sector of the opening. However, 
for calculating stringer stresses in the longitudinal 
interval of the opening, the elementary flexure formula 
appears to be adequate if the bending is symmetrical. 
This is an important result because stresses due to 
symmetrical bending must frequently be superposed 
on stresses due to antisymmetrical bending. Such is 
the case if a fuselage that sustains a vertical tail load 
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Fic. 13. Shearing stresses for Test No. 3 with torque loading. 
Curve for Station 5.5 is derived from Eq. (20). 
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Fic. 14. Shearing stresses for Test Nos. 4 and 7 with torque 
loading. Curve for Station 5.5 is derived from Eq. (20). 


has a rectangular opening that is not centered on a 
horizontal or vertical radius of the cross section. 


Example 


To illustrate the application of the results of the fore- 
going investigation, a typical stress-analysis problem 
is presented. An airplane is subjected to a tail load 
which creates torsion and shear on cross sections in 
the region of a door in the port side of the fuselage 
(Fig. 19). The skin thickness is 0.025 in. The 
stringers are zee-sections having 0.0676-sq.in. cross- 
section area and spacing of 5.6 in. Hence: 


90000 IW. LB 
© TEST NOS 
4 TEST NO.6 
~ 
a: 
4 M/2At 
2000 8 2 
@| | s7a. 
° 
© 
STA| SS | | | 
5 
M/2At 
-/20 -60 60 420 
@ DEGREES 
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loading. Curve for Station 5.5 is derived from Eq. (20). 
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Fic. 16. Shearing stresses for Test Nos. 8 and 9 with torque 
loading. Curve for Stations 4.5 and 5.5 is derived from Eq. (20). 


t. = 0.025 + (0.0676/5.6) = 0.0371 in. 


The opening includes an angle of 62°, and therefore 
B = 149° = 2.60 rad. The radius of the fuselage at 
the mid-section of the cutout is r = 53.4 in. The 
cross-section area of either longeron is a = 0.35 sq.in. 
The length of the opening is 22 = 33 in. 

As in Fig. 1, the positive x-direction is taken to be in 
the direction of the opening. The ordinate y and the 
angle 6 must then be measured as shown in Section A-A, 
Fig. 19, in order to conform with Fig. 1. The agree- 
ment between these figures may be seen by the fact 
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that Section A-A, Fig. 19, may be rotated into co- Fic. 18. Stringer stresses for cantilever test in which the 


It is seen that according to the opening was on the tension side of the cylinder. All curves are 


incidence with Fig. 1. derived from the elementary beam formula, ¢ = My/I. 


sign conventions of Fig. 1, the indicated shear is 
positive and the torque is negative—i.e., V = 30,000 
Ibs., M = —10°in.lbs. Likewise, the bending moment 
at the mid-section of the opening, as determined by the 
30,000-1b. down-tail load, is negative, since it tends to 
rotate the section in the way that a right-hand screw 
would rotate in regressing on the x-axis. The bending 
moment at the mid-section A-A, which is the moment 
of the 30,000-Ib. tail load about that section, is specified A 
to be M,, = —2,250,000 in.lbs. 

: The essential constants are now evaluated. 


h = 0.35/(53.4 X 0.0371) = 0.176 


k = V6 X 53.4? X 0.0871 X 1.3/16.52 X 0.025 = 11.0 


— 2,250,000 - 
"e 53.4? (2.60 + 0.515 X 0.857 + 2 X0.176 X 0.515?) 
= —252 
2(1 + k?)-(sin 8B + k cos B) = —0.146 
B+ sin Bcos B = 2.16 
2sin 8 = 1.03, 1/k = 0.091, 1+ Ak = 2.94 It follows 

cos — Asin B = —0.948 
q oa V/r = —562, M/2r? = —175.4. T= g/t = —3,160e — 10,400 cos 6 — 540 (A) 
o = N/t, = 438ze!!® © + (—1312 + 6,800) sin 0 (B) 


At the forward and aft ends of the cutout, z takes 


30000 LB 


Fic. 19. Diagram of airplane used for illustrative example. 


With these numerical values, Eq. (19) becomes 
—0.146A; + 2.1642 + 1.03A3; = —562 


0.0914; + 0.515As + 2.604; = —175.4 the respective values —16.5 in. and +16.5in. Hence, 
2.944, — 0.94842 + A; = 0 Cae = 7,240e°°-” + 4,640 sin 0 (C) 

The solutions of these equations are Stwa. = —7,2400 + 8,960 sin 6 
A; = —79.0, A, = —260, A; = —13.4 For negative values of 6, the signs of the first terms 
; of Eq. (C) are reversed, as is indicated by the factor 

Then, for @ > 0, Eqs. (20) and (21) yield (1) in Eq. (21). 

q = —79e¢ - ®) _ 260 cos @ — 13.4 Stress curves may be plotted from Eqs. (A) and (C). 


N = 16.26ze° — © + (—4.862 + 252) sin 0 Standard tables of the function e~* determine values 
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TABLE 1 
6 11 (6 — 8) (Rad.) e11(@ — B) —3,160e11(@ — 8) — 10,400 cos 6 t (by Eq. (A)) 
149° —28.6 0 0 — 10,400 — 10,940 
15° 134° —25.7 0 0 — 10,040 — 10,580 
30° 119° —22.8 0 0 —9,000 —9,540 
45° 104° —20.0 0 0 —7,360 —7,900 
60° 89° —17.1 0 0 —5,200 —5,740 
75° 74° —14.2 0 0 —2,700 —3,240 
e 59° —11.3 0 0 0 —540 
105° 44° — 8.4 0 0 2,700 2,160 
120° 4 — 5.6 0.0037 —10 5,200 4,650 
135° 14° — 2.68 0.0684 —220 7,360 ’ 
140° 9° — 1.73 0.178 —560 7,960 6,860 
145° 4° — 0.77 0.463 — 1,460 8,520 6,520 
149° 0 0 1.000 —3,160 8,900 6 
of the factor e@°-”. Tabular values are obtained | 50000 iB 
from Eq. (A) and are shown in Table 1. fo" tM, LB. 
The axial stresses o may be similarly tabulated by Myy=-2.25 X 10° NLA 
means of Eq. (C). The calculated stress curves are x 3 
shown in Fig. 2. Although a symmetrical bending «| SHEARING STR, 
moment due to a load on the fin and rudder was not 9000 ~~ 
included in this example, such a condition could be pon ak » 
incorporated in the calculations by superposing an o 
ordinary flexural stress distribution (ratio of bending r 
moment to section modulus) on the stringer stresses N- A. sramesr |_| 
wa "4 STRESS 
that are plotted in Fig. 20. ~ 8000 | I iL 
STATION AT ‘AFT 
4000 END OF CUTOUT, 
CONCLUSIONS 
The theory has been restricted to the consideration 
of a segment of an open beam that is in equilibrium RN 
under the action of forces and moments that are ap- S > 
plied at the ends. The traction on either end is such 
that the resultant force is perpendicular to the axis of STATION\ AT 
= FORWARD END 
symmetry of the cross section, and such that the re- OF CuTouT, 
sultant moment vector lies in the plane of symmetry ee 


of the beam. With the assumption that hoop stresses 
in the skin are negligible, it has then been demonstrated 
that minimum strain energy results when the stresses 
throughout the beam are determined by Eqs. (19), 
(20), and (21). In particular, Eqs. (20) and (21) show 
that the shearing stress does not vary with the axial 
coordinate z, and that the axial stresses at the mid- 
section (g = 0) are determined by the elementary 
flexure formula (ratio of bending moment to section 
modulus). 

Tests that were conducted for the purpose of esti- 
mating the applicability of these results to the stress 
analyses of cylinders with cutouts have furnished the 
following conclusions: 

(1) In the cutout interval the shearing stress does 
not vary appreciably with the axial coordinate z. 
(See Figs. 16 and 17.) 

(2) Pure torsion creates no axial stress at the mid- 
section of the cutout interval. 

(3) In the case of lateral bending which produces 
a symmetrical axial-stress distribution about the x- 
axis, the elementary flexure formula (ratio of bending 
moment to section modulus) applies accurately for the 
cutout interval. (See Fig. 18.) 


Fic. 20. Stress curves for portion of fuselage in interval of 
the cargo door (see Fig. 19). Curves are derived from Eqs. 
(20) and (21). 


(4) In the case of bending which produces an 
antisymmetrical stress distribution about the x-axis, 
the elementary flexure formula determines the axial 
stresses at the mid-section (z = 0) of the cutout interval. 
(See Figs. 9 and 10.) 


(5) In the case of bending which produces an anti- 
symmetrical stress distribution about the x-axis, the 
stresses in the region opposite to the cutout are deter- 
mined accurately by the elementary flexure formula. 
However, this formula is inapplicable for large values 
of the angular coordinate 0; i.e., for values of |6| 
approaching 8. (See Fig. 6.) 


(6) When the beam cross section is symmetrical 
about the x-axis, the axial stresses and the shearing 
stresses due to combined torsion and antisymmetrical 
shear and bending are respectively odd and even func- 
tions of the angular coordinate 6. 
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(7) The stress distribution varies rapidly in the 
cutout interval as the longerons are approached. 
High axial stresses usually occur at the longerons near 
the corners of the opening. (See Figs. 5 to 17.) 

(8) Eqs. (20) and (21) are fairly accurate for de- 
termining the stresses in the cutout interval of a sym- 
metrical cylindrical beam of constant skin thickness 
and constant equivalent plate-stringer thickness. (See 
Figs. 5 to 17.) Although the case of variable thickness 
was not investigated, a reasonable approximation pro- 
cedure is to calculate N and q with average thick- 


1947 


nesses. The stresses may then be calculated by the 
formulas ¢ = N/t, and tr = q/t, in which ¢ and f, are 
local values of the skin thickness and the equivalent 
thickness. 

(9) Some persistence of the stress distribution in 
the region of the opening is observed in the bays ad- 
jacent to the opening. (See Figs. 5 to 18.) However, 
the data indicate that the transition to the undistorted 
stress pattern is fairly rapid, and therefore hoop 
stresses must be significant in the bays beyond the 
cutout. 


Appendix 


It has been remarked that, when a segment of an open cylindrical beam is loaded by self-equilibrating forces 
that are applied at the ends, the minimizing of the strain energy furnishes the result that the axial stresses at the 


mid-section (zg = 0) are distributed according to the elementary flexure formula. 
of this statement for the case in which the beam loading is antisymmetrical. 


This appendix supplies a proof 
The case of symmetrical loading is 


covered by the elementary beam theory that is given in texts on mechanics of materials. 
In order to take account of the auxiliary conditions, Eq. (6), the Lagrange multiplier method is again employed, 
The strain energy expression, Eq. (7), is then modified by the addition of a linear combination of the integrals of 


Eq. (6); i.e., 


0 No + 2(1 + v) (t.q?/t) + oiNo + sin + o3No cos 6\dé + 


The Euler equation, in this case, reduces to the 
statement that the derivative of the integrand with 
respect to Np is zero. Consequently, with a slight 
modification of the ¢’s there results 


No = o1 + oe sin 6 + 3 cos 0 (A) 
When this expression is substituted in Eq. (6), the 
integrals may be evaluated, and the following equations 
are obtained: 
+ $3 sin B = 
$28 — $2 sin B cos B = (—Myo/r?) — 2(To/r) sin B 
2¢; sin B + $38 + ¢3 sin B cos B = 0 
It follows 
oi = $3 = 0 
— sin B 
r?(8 — sin 8 cos B) 


= 


Consequently, Eq. (A) yields 


(Myo + sin B) sin 6 


No = 
r?(8 — sin B cos B) 


(B) 


Now, by Eq. (2), 
To = (C) 
Consequently, Eq. (B) yields 


+ sin 8) sin B 


To = 
r(8 — sin B cos B) 


Solving for 7» gives 


sin B 


To 
r(8 — sin B cos B + 2) sin? B) 


Hence, by Eq. (C) 


— M, sin B 
r?(B — sin B cos B + 2) sin? B) 


N,|° = (D) 


Now, by Eq. (B), the value of Ny is proportional to 
sin 6; i.e., the following proportion exists: 

No:sin = N)|*:sin B 
Therefore, Eq. (D) yields 


— M,» sin 0 
(E) 
r?(8 — sin B cos B + 2d sin? B) 
This is identical to Eqs. (9) and (10), which are a mathe- 
matical expression of the ratio of the bending moment 
to the section modulus. 
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The Use of Damping Screens for the 
Reduction of Wind-Tunnel Turbulence 


HUGH L. DRYDEN* anp G. B. SCHUBAUERTt 
National Bureau of Standards 


SUMMARY 


The effect of a damping screen in reducing the turbulence of 
an air stream is found to depend on the pressure-drop coefficient, 
k, of the screen. According to a plausible theory, which has 
been experimentally confirmed, the intensity of the turbulence is 


reduced in the ratio 1/VV 1 + & by a single screen and in the ratio 
1/(1 + R)"/2 by x identical screens in series. These relations 
show that for the same power expenditure a greater turbulence 
reduction can be obtained with a number of screens of small k 
than with a single screen of coefficient mk. Factors that place a 
lower limit on velocity fluctuations are turbulence produced by 
the screens themselves and tunnel noise. 


INTRODUCTION 


i ORDER THAT THE AIR STREAM in a wind tunnel may 

act on a body in the same way that it would if the 
body were moved through still air, it is necessary that 
the stream be uniform and steady in speed and direc- 
tion. This ideal can never be perfectly attained. 
There are always present small eddies of varying size 
and intensity which are collectively described as the 
turbulence of the air stream. These eddies arise from 
the wind-tunnel fan, from the guide vanes at the corners 
and from the wind-tunnel walls upstream from the 
working section. In a satisfactory wind tunnel the 
turbulence conditions at any location are statistically 
constant, and it is possible to define the state of the 
turbulence in quantitative terms. The average inten- 
sity of the small eddies can be stated in terms of the 
root-mean-square values of the components of the ve- 
locity fluctuation. An average size of the eddies can be 
determined from measurements of the correlation of 
the velocities at neighboring points as a function of 
their distance apart. 

Even though the turbulent fluctuations amount to 
less than 1 per cent of the mean speed, many aerody- 
namic measurements are greatly dependent on the 
exact values of the intensity and scale of the turbulence 
present. In free flight it has been found that there are 
no disturbances of sufficiently small scale to produce 
appreciable aerodynamic effects—i.e., the turbulence 
of the atmosphere may be regarded as zero. Hence, 
the designers of wind tunnels strive to reduce the in- 
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tensity of the turbulence in wind tunnels as much as 
possible. 

The design and construction of wind tunnels of low 
turbulence is a comparatively modern development for 
two reasons. First, the effect of increased turbulence 
in many aerodynamic measurements is qualitatively 
the same as an increase of Reynolds Number (i.e., in 
size of model or air speed). For a time, an increased 
turbulence was considered desirable as giving an in- 
creased ‘“‘effective’’ Reynolds Number. Second, until 
recently no effective method was known of reducing 
the wind-tunnel turbulence to the point where there 
was no perceptible effect of turbulence on the measure- 
ments. About 10 years ago the intensity of the turbu- 
lence in the average wind tunnel was in the range 0.5 to 
1.0 per cent of the mean speed. Today, there are 
several wind tunnels in which the intensity is of the 
order of 0.02 to 0.05 per cent. 

During the period about 15 years ago, when in- 
creased turbulence was considered to be possibly ad- 
vantageous, screens of coarse mesh were used for in- 
creasing the turbulence of the air stream. It may seem 
paradoxical to the uninitiated that now screens are 
used for decreasing the turbulence. At this point an 
analogy may be helpful, for the ‘“‘smoothing”’ of an air 
stream has some features in common with the ‘‘smooth- 
ing’ of a metal surface. Thus, to smooth a metal 
casting the machinist may begin with coarse file or 
rasp that takes out the gross irregularities and intro- 
duces a new roughness of finer scale. This process is 
repeated with finer files, coarse and fine emery cloth, 
etc., until the desired smoothness is obtained. In 
somewhat similar fashion a screen ‘“‘smooths’’ an air 
stream to the extent of decreasing turbulent motions 
of larger scale than the mesh size while at the same time 
introducing turbulent motions of smaller scale. For- 
tunately, the small-scale turbulence decays much more 
rapidly than the large-scale turbulence, and, hence, at 
some distance from the screen the overall effect is a 
smoothing of the airflow in which both the intensity 
and the scale of the turbulence are reduced. Damping 
screens differ from turbulence-producing screens by 
being of finer mesh and by being located in the large 
cross section ahead of the entrance cone at a great dis- 
tance from the working section rather than in the 
working section itself. 

In the course of the investigations of turbulence and 
boundary-layer flow, carried out at the National Bureau 
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of Standards in 1940 with the cooperation and financial 
assistance of the National Advisory Committee far 
Aeronautics, considerable information was obtained on 
the use of damping screens for the reduction of wind- 
tunnel turbulence. A theory to explain the action of 
damping screens was formulated, and theory and ex- 
perimental results were summarized in an N.B.S. Re- 
port dated August 1, 1940, to the N.A.C.A. Publica- 
tion of the results has been delayed because of the need 
for wartime security; but even though the delay has 
been considerable, the work of 1940 is the only known 
systematic investigation of damping screens in which 
hot-wire measurements of turbulence were made. The 
results may therefore be of interest, since damping 
screens are widely used in wind tunnels of low turbu- 
lence. 

The use of screens in wind tunnels to improve the 
flow dates from some of the earliest wind tunnels. 
Prandtl! called attention to the use of screens for im- 
proving the velocity distribution and gave a simple 
theory of the effect. The use of screens to reduce tur- 
bulence is, as already noted, of fairly recent origin. 
In 1934 a reduction in turbulence was noted behind a 
screen of 8 meshes per in. in the old N.B.S. 4!/.-ft. tun- 
nel. In the N.A.C.A. smoke tunnel at Langley Field, 
cloth was used over the tunnel entrance to secure 
steady flow. Nevertheless, the possibilities of damping 
screens were not fully realized until 1938 when a 60- 
mesh screen was installed in the entrance cone of the 
old N.B.S. 41/2-ft. tunnel and turbulence measurements 
made in the test section showed a reduction in turbu- 
lence from 0.54 to 0.19 per cent. 

Various theories of the effect of screens have been 
proposed in addition to Prandtl’s, among which are 
Collar’s* regarding the effect of a screen on the velocity 
distribution and Batchelor’s* regarding the effect of 
screens on both velocity distribution and turbulence. 
The relatively simple theory given herein is the only 
known one that agrees even approximately with the 
observed effect of screens on turbulence. Although the 
experiments antedated the successful theory, the theory 
is discussed first. 


SYMBOLS 


A = cross-sectional area of wind stream, or duct 

U = mean speed at any point 

p = mean static pressure at any point : 

“u = component of velocity fluctuations produced by turbu- 
lence, parallel to mean flow, measured with respect to 
mean speed 

v,w = mutually perpendicular components of velocity fluctua- 
tions produced by turbulence, normal to mean flow 

u’, v’, w’ = root-mean-square values of u, v, w 

U’ = mean turbulence intensity defined by U’ = 
V 1/s(u't + + w’?) 

p = density of the air 

E’ = kinetic energy of turbulence = !/29(u’? + v’? + w’?) = 

k = pressure-drop coefficient of a screen 

n = number of screens 
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d = wire diameter of a damping screen or rod diameter of a 
turbulence-producing grid 

L = scale of turbulence 

x = distance downstream from a screen 

t = time fora particle to travel from a screen to some point x 

C =contraction ratio of a wind tunnel, A,/A, (see subscripts) 


Subscripts 


t = test section of a wind tunnel 

s settling chamber of a wind tunnel 

1,2,3 = sections upstream, at a screen, and downstream 
from a screen, respectively, in a duct of con- 
stant cross section (see Fig. 1) 

conditions at a particular time—used in decay 
formulas 


THEORY OF EFFECT OF SCREENS ON TURBULENCE 


Consider a straight duct of constant cross section 
with a screen placed entirely across section 2 as shown in 
Fig. 1. Assume now that a turbulent air stream is 
flowing in the duct with constant mean speed U. One 
may regard the turbulence as motions of the air in all 
directions, with component velocities u, v, and w with 
respect to a coordinate system moving with the mean 
speed U. The mean kinetic energy per unit volume of 
the turbulence is then 


= 0? + (1) 
The energy per unit volume of the mean motion is 
b + '/2p UV? 


Let us consider for the present that the screen is so 
fine that it may be regarded simply as a plane of uni- 
formly distributed resistance across which some of the 
energy of the flow is converted into heat. It is assumed 
that the speed is well down in the subsonic region where 
temperature changes across the screen are so small that 
they have a negligible effect on the properties of the air. 
The problem now is to find the dissipation of energy 
per unit volume across the screen. 

The resistance of a screen is proportional to the 
density of the air and the square of the speed. The 
resistance per unit area may be written k('/2)pV’, 
where k is commonly called the pressure-drop coeffi- 
cient and V is the speed close to the screen but not 
within the pores. It is important to bear in mind that 
V may be different from the speed at a great distance 
from the screen if changes occur in the cross section of 
the stream as a whole or in the cross section of any 
elementary stream tube. Since the mean speed U is 
uniform and the duct has constant cross section, U 
cannot change anywhere in the duct. Hence, the re- 
sistance per unit area of screen to the mean flow is 


R(*/2)pU? 


and this is also the energy loss per unit volume. 
fore, in the terminology of Fig. 1 


R(*/2)pU? = pi + — ps — 


and since U; = U3 
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k'/2pU? = pi — ps (2) 


Since the turbulent velocities are correlated only 
over a small portion of the duct area, it is not required 
that E’ remain constant. In this case changes may 
occur in the cross section of the stream tubes with which 
E' is associated. We therefore write the mean resist- 
ance per unit area of screen to the turbulence as kE,’ 
assuming that the value of k is the same as for the mean 
flow. 

Since this is also the energy loss per unit volume 


kE,’ = Ey’ E;’ (3) 


The absence of pressure terms in Eq. (3) is the result of 
the concept that the turbulence pattern is carried along 
with the mean flow but is not forced through a resist- 
ance by it. In other words, there is no transfer of 
energy to the turbulence from the mean flow. 

It will now be assumed that the changes in E’ occur 
on the approach side of the screen and that the E’ at 
the screen is different from £,’. Since all of the dis- 
sipation is assumed to take place in the plane of the 
resistance, no further reduction of energy will occur 
after the flow has passed the resistance. Therefore 
E,' = E;’ and Eq. (3) reduces to 


E,'/E;' = 1/(1 + k) (4) 


Eq. (4) expresses the reduction of turbulent energy 
by a screen of coefficient k. As mentioned previously, 
k is usually termed a pressure-drop coefficient because 
of its role in Eq. (2). Eq. (2) is in fact the relation used 
to determine k from the measured pressure drop ~; — 
p;and the measured speed U. 

By experiment the root-mean-square components 
of the turbulence wu’, v’, and w’ are determined sepa- 
rately, and E’ is given by '/2 p times the sum of their 
squares. It is convenient, however, to use the sum of 
the squares to derive a mean turbulence intensity U’, 
defined by 


so that if the turbulence is isotropic u’ = v’ = w’ = U’. 
The relative intensity is then expressed as U’/U, where 
Vis any convenient reference speed. 

It will be seen that Eq. (4) applies also to U’?, and 
the fractional reduction in U’ produced by a screen is 


= 1/Vi+k (5) 
Ko 
' P, 
| 
E, Ee 
Us 


Fic. 1. Schematic diagram of a damping screen in a duct of 
constant cross section. 


The simplifications and assumption involved in the 
foregoing development furnish ample grounds for 
questioning the validity of the answer. About all 
that can be said in defense of the methods is that the 
process can take place as assumed and that the assump- 
tion of different processes would generally lead to ab- 
surd results. For example, if EZ,’ is set equal to E,’ or 
to some value intermediate between EF,’ and E;’, nega- 
tive values are obtained for E;’/E,’ for sufficiently large 
values of k. 

The assumption that the energy lost at the screen 
goes directly into heat energy is known to be incorrect. 
The energy goes first into small jets through each pore, 
then into tiny eddies and heat, and finally into heat as 
the eddies decay. The eddies constitute high turbu- 
lence, and just back of a screen the turbulence may be 
extremely high, much higher in fact than it would 
generally be without a screen. However, if the screen 
is fine, the scale of the turbulence is small and the eddies 
decay rapidly. 

A considerable amount of experimental and theoreti- 
cal information exists‘ on the scale Z and the inten- 
sity u’ of the turbulence produced by square-mesh 
grids. The turbulence produced by grids is found to 
be isotropic, so that u’ = v’ = w’ = U’. When 
lengths are expressed in terms of the rod diameter d of a 
grid, it is found that L/d and U’/U have approximately 
the same values for all grids at the same number of rod 
diameters downstream in a stream of constant mean 
speed. At 212.5 rod diameters, for example, L/d is 
about 1.5 and U’/U is about 0.02, This distance 
amounts to 1.5 in. for a screen with wires 0.0075 in. in 
diameter. It is only at distances of several feet that a 
net beneficial effect would be expected. 

The experimental results for the decay of turbulence 
behind grids of wires 0.05 to 0.2 in. in diameter can be 
represented approximately by the relations 


(Uo’/U')? = 1 + 0.58 [Uo'(t—to)/Lo] (6) 
(L/L)? = 1 + 0.58 [Uo'(t—to)/Lo] (7) 


where Up’ and Lp» are the values of the mean turbulent 
intensity and the scale at a time f after a particle travel- 
ing with the mean motion passes the grid and where 
U’ and L are the values at time ¢. It is not known how 
accurately these relations represent the decay behind 
screens of wire diameters 0.007 to 0.02 in. They prob- 
ably give conservatively high values of U’, for it is 
known that as the wire diameter decreases a Reynolds 
Number is eventually reached at which no eddies are 
shed from the wires. Thus a screen made up of suf- 
ficiently fine wires should create no turbulence of its 
own. 

Returning to the consideration of Eq. (5), it may be 
noted that Prandtl,’ in dealing with variations in mean 
speed over the cross section of a stream, obtained re- 
ductions in the differences from point to point in the 
ratio 1/(1 + k), which is the square of the value found 
for the reduction of turbulence. The difference in the 
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expressions arises from the fact that the average values 
of the velocity variations over a long time are zero for 
turbulence but not for differences of mean speed across 
the section. 

The price paid for the reduction of turbulence by a 
damping screen is a loss of energy equal to the pressure 
drop ~i — p; times the total volume of air moved 
through the screen. The pressure drop through a 
single screen is k('/2)pU,”, where U, is the speed at the 
screen. In a wind tunnel the speed at the screen can 
be diminished by locating the screen at the largest 
cross section of the entrance section. The speed can 
be reduced at will by increasing the contraction ratio; 
it is inversely proportional to the contraction ratio. 

It is obviously advantageous to use a number of 
screens, say 7, of coefficient & in series rather than a 
single screen of coefficient nk, both giving the same en- 
ergy loss, since in the first case the reduction of tur- 
bulence is by a factor 1/(1+k)"/” as compared with 1 + 
(1+nk)'” in the second. Thus, if n is 5 and k is 2, the 
factor with five screens is 0.064 as compared to 0.301 
for a single screen of the same total resistance. 

The theoretical ratios of the turbulence behind the 
screens to that ahead of the screens for several values of 
n and k are given in Table 1. 

It should be repeated that these relations take no 
account of turbulence set up by the screen itself, and 
they cannot be expected to hold unless the turbulence 
produced by the screen has decayed to insignificant 
proportions at the point where an observation is made. 
The lower limit of the turbulence level attainable by 


TABLE 1 


Theoretical Reduction of Turbulence by » Damping Screens of 
Pressure-Drop Coefficient k 
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the use of damping screens is determined by the fine- 
ness and uniformity of the screen. Departures from 
theory are therefore to be expected when the turbu- 
lence is low and when large reductions are involved. 
There is also another factor affecting comparison with 
theory which becomes increasingly important as the 
turbulence is reduced—namely, the apparent turbu- 
lence due to sound in a wind tunnel. 


EXPERIMENTAL RESULTS 


The investigation was conducted in the rebuilt N.B.S, 
4'/,-foot wind tunnel shown in Fig. 2. The damping 
screens were installed in the settling chamber, not 
necessarily in the positions indicated in tlie figure, and 
turbulence measurements were made with the hot-wire 
turbulence-measuring equipment described in reference 
5. Briefly, the method employs a single platinum 
wire 0.0002 in. in diameter about 0.2 in. long normal 
to the wind for the measurement of u’ and a pair of 
such wires in the form of an X for the measurement of 
v’ and w’. The plane of the X was horizontal for the 
measurement of the horizontal component v’ and ver- 
tical for the measurement of the vertical component w’, 
Fluctuating voltages from the wires heated by a con- 
stant current were amplified and compensated for ther- 
mal lag, and the mean-square was then measured by a 
thermocouple-type milliammeter. Using calibration 
data for wires and amplifier, values of u’, v’, and w’ 
were calculated. . The bulk of the measurements were 
made at a single station in the test section about 6 ft. 
from the upstream end. Some measurements were 
also made in the settling chamber. 

Data are available on four single screens, two screens 
in series, three screens in series, and six screens in series. 


nk = 5 nk = 10 nk = 15 nk = 20 Neither the position in the settling chamber nor the 
n Ratio = 1/(1 + k)»/? spacing between screens had a measurable effect on the 
A p 
results. The spacing was usually between 6 and 12 in., 
3 0.230 0.111 0.068 0.047 but to study the effect it was varied in one case from 2 to 
; 28 in. The detailed results of the measurements in the 
6 0.162 0.052 0.0234 0.0123 test section are given in Table 2 and of the measurement 
ph bens in the settling chamber in Table 3. To facilitate com- 
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Fic. 2. Elevation view of 4'/,-ft. tunnel. 
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TABLE 2 
Turbulence of N.B.S. 41/.-Ft. Wind Tunnel with Various Damping Screens 


Measurements at a Single Station in the Test Section) 


I. No Damping Screens 


U;, % U1, U1, 
Ft. per Sec U,’ Ft. per Sec. ia Ft. per Sec a” 
31.6 0.335 32.4 0.373 
45.9 0.090 44.8 0.333 45.7 0.363 
67.1 0.094 65.2 0.302 66.4 0.364 
87.4 0.093 87.6 0.289 87.6 0.326 
108.8 0.089 109.6 0.301 109.7 0.335 
133.0 0.084 131.7 0.317 132.5 0.346 
Us, Us Us, Us Us, w’ 
Ft. per Sec. Ft. per Sec. Ft. per Sec. ii Ft. per Sec. Be. Ft. per Sec U, “ 
II. One 18-Mesh Screen, 0.011-In. Wire 
9.41 62.1 .050 we es 31.8 0.195 30.7 0.217 
15.20 100.3 .917 44.9 0.041 44.7 0.195 47.2 0.216 
20.0 132.0 . 842 64.9 0.046 65.4 0.210 63.1 0.226 
87.0 0.048 88.5 0.201 88.3 0.208 
108.6 0.052 109.2 0.199 110.7 0.239 
136.1 0.054 136.5 0.198 132.8 0.235 
III. One 20-Mesh Screen, 0.017-1n. Wire 
8.82 58.2 2.572 29.5 0.026 28.3 0.185 28.8 0.187 
15.35 101.3 2.220 42.0 0.031 41.5 0.171 44.4 0.183 
19.60 129.3 2.069 58.6 0.033 57.2 0.177 60.3 0.186 
81.6 0.039 83.0 0.162 80.8 0.185 
104.6 0.041 104.7 0.172 104.8 0.199 
127.0 0.051 127.0 0.163 125.7 0.199 
IV. One 24-Mesh Screen, 0.0075-In. Wire 
9.54 63.0 0.904 31.1 0.039 30.9 0.176 31.4 0.200 
14.93 98.6 0.747 40.9 0.039 40.5 0.189 43.8 0.224 
19.52 128.8 0.705 60.4 0.045 60.0 0.192 60.4 0.206 
83.4 0.047 84.0 0.194 83.6 0.219 
104.8 0.053 105.3 0.228 105.3 0.244 
130.5 0.058 131.1 0.245 132.8 0.259 
V. One 60-Mesh Screen, 0.007-In. Wire 
9.50 62.7 5.300 30.2 0.019 28.9 0.075 29.5 0.083 
14.90 98.4 4.356 46.1 0.025 44.7 0.108 45.7 0.109 
19.70 130.0 .903 62.6 0.029 62.8 0.122 62.8 0.116 
88.4 0.033 89.0 0.134 88.0 0.129 
109.4 0.038 110.6 0.141 108.8 0.133 
Ut %, Ut, % Us, w 
Ft. per Sec 1 i Ft. per Sec. U?’ Ft. per Sec. ug * 
Two 18-Mesh Screens, 0.011-In. Wire. Spacing 6 In. 
30.7 0.022 28.9 0.117 30.3 0.125 
46.8 0.026 45.1 0.125 46.2 0.130 
63.2 0.031 61.4 0.136 63.2 0.147 
86.4 0.034 87.1 0.144 86.6 0.155 
108.6 0.041 109.0 0.166 110.0 0.155 
132.4 0.190 132.0 0.157 
Three 18-Mesh Screens, 0.011-In. Wire. Spacing 6 In. 
30.0 0.016 30.4 0.071 28.1 0.077 
42.3 0.020 45.0 0.092 43.5 0.087 
62.4 0.026 61.1 0.112 61.0 0.099 
87.6 0.032 87.2 0.111 87.8 0.118 
108.2 0.037 109.1 0.114 109.7 0.123 
131.2 0.043 129.3 0.119 130.0 0.121 


(Three 20-Mesh, 0.017-In. Wire; Three 24-Mesh, 0.0075-In. Wire. 
20, 24, 24, 20, 20, 24. 


Spacing 4, 6, 8, 12, 23 In.) 


Order of Screens, Upstream 


to Downstream, 


28.8 0.019 27.5 0.009 28.5 0.011 
41.7 0.021 40.4 0.020 41.9 0.020 
60.2 0.027 60.6 0.026 61.6 0.021 
84.0 0.034 84.5 0.037 86.2 0.036 
107.4 0.039 108.1 0.044 109.2 0.045 
123.1 0.046 124.0 0.047 124.0 0.048 


parison with theory the results in the test section per- 
taining to a single speed (about 109 ft. per sec.) have 
been summarized in Table 4. The values shown in 
the seventh column of Table 4 agree with the predic- 
tions of the theory well within the precision of the ex- 
periments, except for the six-screen combination. It 


will be shown later that the lack of agreement for the 
six-screen combination may be attributed partly to 
incomplete decay of the turbulence produced by the 
last screen and partly to an effect of sound, discussed 
in the next section. It is known from a few traverses 
that there were small variations across most of the 
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TABLE 3 
Turbulence in Settling Chamber of N.B.S. 4!/.-Ft. Wind Tunnel 


(Measurements at a Single Station. 


Experimental Value of U;/U, = 6.6) 


Ui, ul U1, U1, w’ w’ 
Ft. per Sec. U,’ % U? % Ft. per Sec. JU,’ % U? % Ft. per Sec. JU,’ % U? % 
I. No Damping Screens 
33.0 1.51 0.229 31.9 1.57 0.238 33.4 1.67 0.253 
44.6 1.69 0.256 45.5 1.58 0.239 45.5 1.62 0.245 
62.3 1.69 0.256 66.8 1.47 0.223 64.6 1.58 0.239 
84.1 1.94 0.294 86.6 1.49 0.226 85.8 1.52 0.230 
105.0 1.97 0.299 107.2 1.51 0.229 107.2 1.56 0.236 
127.9 2.31 0.350 131.3 1.47 0.223 131.3 1.65 0.250 
II. Six Damping Screens. (Arrangement VIII in Table 2. Measurements 251/2 In. from Downstream 24-Mesh Screen) 
28.9 0.051 0.008 29.1 0.036 0.005 29.3 0.037 0.006 
39.1 0.090 0.014 39.7 0.040 0.006 39.7 0.041 0.006 
58.2 0.139 0.021 58.2 0.126 0.019 58.2 0.130 0.020 
81.1 0.148 0.022 80.0 0.131 0.020 79.7 0.180 0.027 
103.7 0.175 0.027 104.0 0.146 0.022 104.0 0.222 0.034 
119.3 0.207 0.031 118.7 0.144 0.022 118.4 0.247 0.037 
39.4 0.054 0.008 
51.2 0.108 0.016 
TABLE 4 


Comparison of Observed and Predicted Reduction of Turbulence by the Use of Damping Screens in N.B.S. 41/.-Ft. Wind Tunnel 


Screen nk 1/(1 + k)n/2 

None 

18-mesh, 0.011-in. wire 0.895 0.726 
20-mesh, 0.017-in. wire 2.18 0.561 
24-mesh, 0.0075-in. wire 0.730 0.760 
60-mesh, 0.007-in. wire 4.20 0.438 
Two 18-mesh, 0.011-in. wire 1.790 0.527 
Three 18-mesh, 0.011-in. wire 2.685 0.383 


Six screens (three 20-mesh, 
0.017-in. wire, three 24- 


mesh, 0.0075-in. wire) 7° 0.0776t 


U; = 109 ft. per sec. 


U’ 

u’ v’ w’ Observed 
0.089 0.301 0.335 0.265 

0.052 0.199 0.239 0.182 0.192 
0.041 0.163 0.199 0.150 0.149 
0.053 0.228 0.244 0.195 0.201 
0.038 0.141 0.133 0.114 0.116 
0.041 0.166 0.155 0.133 0.140 
0.037 0.114 0.123 0.099 0.102 
0.039 0.044 0.044 0.043 0.021 


* 3koo + 3k, where ko and ky are the pressure-drop coefficients for the 20-mesh and the 24-mesh screens, respectively. 


+ (1 + 


test section with large increases in the wake of the 
seams of the screens. In a previous test with the same 
six screens, a 20-mesh screen with 0.017-in. diameter 
wires was tried as the final screen. There was con- 
siderable variation across the section in this case, 
which was taken to indicate that a measurable turbu- 
lence was produced by the final screen either because 
the mesh and wire size were too large or because the 
screen was not sufficiently uniform. It is well there- 
fore to consider the decay of turbulence behind the 


. last screen in greater detail. 


As pointed out in the previous discussion, the decay 
formula, Eq. (6), probably predicts high values for the 
turbulence produced by grids of wires as small as those 
used in damping screens. Nevertheless, it is informa- 
tive to calculate the decay by Eq. (6), using Uo’ and 
Ly’ for grids, and thus estimate the maximum turbu- 
lence produced by a damping screen. For grids the 
following values are assumed: 


Uo’ = 0.02* U,, Lo = 1.5d at to = 212.5/U, 


* For different grids Uo’ at this location varied between 0.014 
and 0.030 times U,. 


To obtain U’ it is assumed that the acceleration of 
the stream through the contraction affects only the 
time for decay. The relative intensities may then be 
expressed in terms of the speed in the settling chamber 
U, or the speed in the test section U;. Taking ¢ as the 
time of decay from a point 2 ft. upstream from the 
entrance section to a point 6 ft. into the test section, 
the relative intensities found for the several screens 
are shown in Table 5. 

If the values in the last column of Table 5 are com- 
pared with those in the last column of Table 4, it is 
seen that the turbulence created by the last screen is 
considerably less than the turbulence passed by the 


TABLE 5 


Turbulence in Test Section Produced by Last Damping Screen as 
Computed from Eq. (6) 


Meshes Wire diameter, U! % u' % 
per in. in. U,’ oO Ui’ 0 
18 0.011 0.228 0.035 
20 0.017 0.282 0.043 
24 0.0075 0.188 0.029 
60 0.007 0.182 0.028 
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screen, except for the case of the six-screen combina- 
tion. Here, the above value of U’/U, is comparable 
with the value actually measured. In the case of the 
six-screen combination, if the last 24-mesh screen pro- 
duced 0.029 per cent and the predicted value of the 
turbulence passed by the screen were 0.021 per cent, 
the expected total intensity would be 


V (0.029)? + (0.021)? = 0.036% 


which stands in fair agreement with the 0.043 per cent 
observed. 

The scanty observations made in the settling cham- 
ber demonstrate that Eq. (6) for the decay of turbu- 
lence does not apply too well to damping screens. Un- 
fortunately, observations were not made of the scale 
of the turbulence in the settling chamber and in the 
test section for each combination. However, if the 
same assumptions are made as those used in comput- 
ing Table 5, one finds that the turbulence produced 
by the 24-mmesh screen at a distance of 25'/2 in. is pre- 
dicted to be 0.38 per cent as compared with observed 
values from 0.04 per cent at low speeds to 0.20 per cent 
at high speeds. The Reynolds Number of the screen 
wires ranges from about 35 to 140, which is in the range 
at which eddies first begin to form behind cylinders. 
The Reynolds Number of the turbulence itself is also 
low. At the highest speed the predicted turbulence 
passing through the screen is 0.20 per cent. It ap- 
pears therefore that there is practically no turbu- 
lence produced by the damping screen itself in this 
case. 


EFFECT OF SOUND 


Assuming that the last screen of the six-screen com- 
binations was producing a significant part of the turbu- 
lence measured, a seventh screen made of silk bolting 
cloth was installed as the final screen. Silk bolting 
cloth was used because it was finer than any wire screen 
commercially available. Some additional reduction 
was obtained, but the discrepancy between theory and 
experiment was even greater than before. 

During the course of these experiments it was noted 
that the wave form of the velocity fluctuations did not 
have the randomness usually associated with turbu- 
lence. This led to the conclusion that something in the 
nature of a vibration was being picked up by the hot 
wire. Tests showed that vibration was not being 
transmitted to the wire itself, and the only other con- 
teivable source was sound. Measurements with a 
sound meter connected to a crystal microphone showed 
the sound level to be between 105 and 110 db. at 100 
ft. per sec, depending on where the microphone was 
Placed. Since it was necessary to keep the micro- 
phone out of the wind, the sound level could not be 
measured at the position of the hot wire. Assuming 
plane waves, the root-mean-square particle velocity 
was calculated for 105 and 110 db. above a base level 


of 10-* ergs per cm.* per sec., giving the following 
results: 


Particle Velocity, Ft. per Sec., or Apparent Per- 
centage Turbulence Due to Sound at Wind 


Db. Speed 100 Ft. per Sec. 
105 0.028 
110 0.051 


The foregoing percentages are what might have been 
indicated by the hot wire because of sound alone; 
the actual value indicated, including sound and all 
components of real turbulence, was 0.043 per cent. 
Obviously, the sound measurements were not accurate 
enough to permit even an approximate correction for 
the sound. It is interesting to note, however, that the 
remaining discrepancy [V (0.043)? — (0.036)*] is 
0.024, and this figure is of the same order as those esti- 
mated from sound measurements. 

As a further test to determine whether particle veloc- 
ity of sound was picked up by the hot wire, a frequency 
analysis was made of the output of the hot wire with 
the tunnel running at 80 ft. per sec. with the six-screen 
combination of Table 2—VIII and of the sound in the 
control room under the same conditions. The results 
are shown in Fig. 3. The similarity of the two curves 
leaves little doubt that a considerable part of the hot- 
wire output is derived from sound. 
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Fic. 3. Distribution with frequency of u-component of tur- 
bulence and sound. Wind speed—80 ft. per sec., u’/U; (total) 
= 0.033 per cent. 


According to Batchelor,’ noise in closed wind tun- 
nels has many sources, two important ones being the 
propeller and the so-called drag noise produced by ir- 
regular eddy motion in the boundary layer at the tunnel 
walls. In the present tunnel propeller noise was 
evident in the measurements, and by ear the roar of the 
propeller seemed to be louder than all other sounds. 
In this tunnel therefore, the propeller is evidently 
the major source of sound, and being in the same arm 
as the test section, the sound is readily transmitted to 
the test section. 

It should be pointed out that the screens are still 
partially responsible for the lower limit placed by sound. 
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In the present arrangement, for example, the propel- 
ler had to run faster and faster for a given wind speed 
in order to overcome the increasing pressure drop 
across the screens as more or denser screens were 
added. The net result was that the propeller noise 
increased with each addition of a screen. 


EFFECT OF CONTRACTION 


The application of the theory of the damping effect of 
screens located in the settling chamber to measurements 
made in the test section cannot be expected to be suc- 
cessful unless the effects of decay and contraction in the 
entrance cone are equivalent to multiplication by a 
constant factor. Measurements were made in the 
settling chamber for two conditions, first with no damp- 
ing screens and again with the six-screen combination. 
Values of the turbulent fluctuations are given in Table 
3, both as percentages of the local speed, U;, and as 
percentages of the speed in the test section, U;. At 

/, equal to about 104 ft. per sec., for example, the aver: 
age turbulence in the settling chamber at the point un- 
der consideration, without damping screens, is 0.257 
per cent of U;. That in the test chamber itself is 0.265 
per cent of U,, which is not significantly different from 
that in the settling chamber. With the six-screen 
combination the average turbulence appears to increase 
from 0.028 per cent of U, in the settling chamber to 
0.043 per cent of U, in the test section. 

The ratio between the speed of the core at the test 
section to the average speed in the settling chamber 
was approximately 6.6. The departure of this value 
from the contraction ratio 7.1 represents the effects of 
the boundary layer in the test section. 

These few observations suggest caution in attribut- 
ing any generality to the result that there was little 
effect of contraction on the energy of the turbulence, 
for the time of transit through the contraction is equiva- 
lent to a distance traveled at the settling chamber speed 
of 9 ft. The decay formula predicts a fairly large decay 
(ratio about 0.3 if the initial scale is 0.6 in.), and it there- 
fore is obvious that the effect of the contraction is 
actually to increase the turbulent energy. In the case 
of the 4'/2-ft. wind tunnel the effects happen to balance 
approximately. 

With no screens present the contraction of the air 
stream does exert a selective effect on the three com- 
ponents u’, v’, and w’, u’ being reduced to about one- 
third of its initial value and v’ and w’ increased by 
about 40 per cent. Such selective effects of a contrac- 
tion are discussed by Prandtl! and by Taylor.’ Their 
theories neglect the decay of turbulence and are based 
on regular types of disturbances without statistical 
considerations. Prandtl predicts a decrease of u’ in 
the ratio 1/C and an increase of v’ and w’ in the ratio 


VV 2 where C is the contraction ratio. This result cor- 
responds to an increase of U’ in the ratio 
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Vv (2C/3) + (1/3C?), which equals 2.1 for a contrae- 
tion ratio of 6.6. Taylor obtains a result dependent on 
the assumed form of the disturbance—in some cases a 
decrease in the turbulent energy, in others an increase. 
Neither theory is free from objection. 

The principal advantage of a large contraction ratio 
has already been mentioned—namely, that of permit- 
ting the installation of damping screens without an 
excessive energy loss. Furthermore, if the turbulent 
energy remains constant or decreases, the ratio of the 
mean turbulent intensity U’ to the mean speed de- 
creases because of the increasing mean speed as the 
cross section contracts. 


CONCLUSIONS 


The experiments described in this paper show that 
damping screens materially reduce the turbulence of an 
air stream, and by their aid it is possible to obtain tur- 
bulence levels as low as 0.02 to 0.05 per cent. The 
action of the screens is interpreted in terms of a simple 
theory that assumes that the effect of a damping 
screen is partly to absorb the kinetic energy of the 
turbulence. The theory takes no account of any turbu- 
lence set up by the screen itself and therefore may pre- 
dict values too low. This, of course, depends on the 
fineness and uniformity of the screen. The turbulence 
set up by the screens themselves may be the limiting 
factor in the attainment of low turbulence, but in the 
4'/,-ft. tunnel the principal limiting factor is probably 
the tunnel noise. While velocity fluctuations arising 
from a sound field should not be classed as turbulence, 
they, like turbulence, are generally an undesirable source 
of disturbance and should be reduced as much as pos- 
sible. 

The information available on the effect of a contrac- 
tion on turbulence is not sufficient to permit the exact 
computation of the turbulence levels to be expected in 
wind tunnels of different contraction ratios. 
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The Application of Basic Physiologic Data 
the Design of Respiratory Equipment 


LOREN D. CARLSON* anp VERNER J. WULFFt 
Unwersity of Washington and State University of Towa 


ABSTRACT 


Data on respiratory physiology are assembled and unified into 
general statements that may be of use to engineers in the design 
of respiratory equipment. The varied and apparently inde- 
pendent data concerned with partial pressures of oxygen required 
for adequate arterial blood saturation; pulmonary ventilation; 
respiratory rate; human tolerable resistance to flow; and gas 
exchange in the lungs are collected and presented in a generalized 
framework. 


SYMBOLS 


F = flow of oxygen liters per minute S.T.P.D. 

S.T.P.D. = 760 mm. Hg, 0°C., dry 

B.T.P.S. = body temperature and pressure saturated with water 
vapor at 37°C. 

M.I.F. = maximum instantaneous rate of flow, liters per min. 
B.T.PS. 

J = ventilation or minute volume in liters per min. 
B.T.PS. 

0, = fraction of oxygen required 

P = 760 mm. Hg 

: 2 = ambient pressure, mm. Hg 

7 = 273°K. 

Ty = body temperature, 310°K. 


n = number of breaths per min. 
V, = volume economizer or rebreather in liters 


pO» = partial pressure of oxygen’ in mm. Hg 

pCO, = partial pressure of carbon dioxide in mm. Hg 

R.Q. = respiratory quotient 

‘LF. = instantaneous rate of flow, liters per min. B.T.P.S. 
INTRODUCTION 


DESIGN OF CONTINUOUS or diluter demand 
oxygen equipment for use in combat or passenger 
aircraft requires the establishment and use of certain 
basic physiologic data. Design of more specialized 
equipment, such as pressure demand or pressure suit 
equipment, involves additional pertinent physiologic 
data. At present much of this basic physiologic and 
technical data has been accumulated by research divi- 
sions of the armed forces and under the auspices of the 
O.S.R.D. However, it has not been completely collated 
or organized, nor is it available for general use. The 
use of these data in formulating equipment design has 
necessitated an organization of the data for the specific 
purpose. It is the purpose of this paper to present and 
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elaborate those data which are the bases for the design 
requirements of respiratory equipment. The amount 
of data in this field is large and a review of all work is 
not attempted. Only pertinent quantitative data will 


be reviewed. 


The factors that the engineer must carefully evaluate 
are few in number but the volume of data is large and, 
at times, complexly interrelated. Factors to be con- 
sidered are: 


(1) The pressure of oxygen required to maintain ade- 
quate arterial oxygen. 
(2) Pulmonary ventilation: 
(a) Effect of activity, cold, altitude, flight condi- 
tions. 
(b) Inspiratory instantaneous rates of flow. 
(c) Expiratory instantaneous rates of flow. 
(d) Maximum instantaneous rate of flow. 
(e) Mean instantaneous rate of flow. 
(f) Per cent of total volume drawn in at various 
rates of flow. 
(3) Respiratory rate. 
(4) Resistance to flow tolerable to humans. 
(5) Gas exchange in lungs. 


The data pertaining to these factors may have vary- 
ing significance on the design of the various types of 
oxygen apparatus. Basically, all the data should be 
reviewed for any system. 

Generally, data on any experimental subject are 
presented as an average with sufficient data being col- 
lected to insure significance. However, the average 
figure, if used in the design of respiratory equipment, 
can give an erroneous interpretation of equipment 
value. In no piece of equipment can the design be con- 
fined to the mean or average, for by definition that de- 
velopment excludes practically 50 per cent of the indivi- 
dual users from having adequate protection. A tenable 
thesis is proposed that the requirement for equipment 
should be so designated that 90 per cent of its users will 
be completely protected. Thus, in prediction of the 
duration of an oxygen supply, a value on the basis of 
average figures would give the flier an even chance of 
having sufficient oxygen; on the basis of the thesis just 
stated the chances would be 9 to 1. If the average 
figure and the standard deviation are known, then the 
value to be used is determined as: 


average + (1.3 X standard deviation) (1) 
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In a few cases where necessary a minimum figure may 
be obtained by a corresponding subtraction. 


PARTIAL PRESSURE OF OXYGEN 


The question of the partial pressure of oxygen neces- 
sary to maintain adequate arterial saturation has been 
extensively discussed and adequate basic data pre- 
sented. 2 Calculations to determine the necessary 
fraction of oxygen may be determined by: 


f02 = [pO2 + (pCO2/R.Q.)]/ 
[P, — 47 — pCOz + (pCO:/R.Q.)] (2) 


or by a simpler formula 
fO2 = pO2/(P, — 47) (3) 


The desirability of using the latter formula in equip- 
ment has been ably presented by Boothby and Fergu- 
son.” 

The minimum pO, that is desirable has been indicated 
to be that simulating 5,000-ft. breathing air.'~* This 
then is a minimum value for the equipment; the appara- 
tus must deliver at least this quantity of oxygen at all 
conditions of flow and temperature. Alveolar pO, at 
5,000 ft. is given at 81.6 mm. Hg by Boothby.? The 
standard deviation is 4.5 mm. Hg. To cover 90 per 
cent of the cases in maintaining a pO: (alveolar) at 81.6 
mm. Hg, the pO: for tracheal air in Eq. (3) will be fixed 
at 134 mm. Hg. 


PULMONARY VENTILATION 


Pulmonary ventilation or the amount of gas mixture 
breathed has significance in the design of respiratory 
equipment (mask and regulator) and in the amount of 
oxygen placed in the aircraft. The latter may be dis- 
cussed briefly. The amount of oxygen to be carried in 
aircraft certainly cannot be based on an average nor 
entirely on experimental results. Reports on experi- 
mental flights’ indicate that the average ventilation in 
flight expressed at B.T.P.S. is constant up to 33,000 ft. 
The average value is 14.2 liters per min. B.T.P.S. for 
normal flying conditions. Statistical analysis places the 
90 per cent inclusive figure at 18.4 liters per min. 
Data from actual combat flights® present an average of 
11.4 liters per min. B.T.P.S. and the 90 per cent in- 
clusive figure is 18.2 liters per min. These values are 
for A.A.F. airplanes. For cargo or passenger aircraft, 
the value of 10 liters per min. given by Boothby? may 
be more nearly correct. Addition to this to cover 90 
per cent of the cases would give a value of 12 liters per 
min. B.T.P.S. 

An average active ventilation has also been deter- 
mined experimentally. No combat data are now avail- 
able. The average active figure is presented® as 23.2 
liters per min. B.T.P.S. Ninety per cent of the cases 
are covered by a value of 37.1 liters per min. B.T.P.S. 
From these data it is interesting to note that to cover 
adequately approximately 100 per cent of the cases 
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Fic. 1. Standard inspiratory flow-time curves for 18 and 
82 liters per min. Ventilation reconstructed from data presented 
in this report. Inspiratory flow curves reconstructed from the 
generalized data presented in Figs. 3 and 5. These inspiratory 
curves will result in minute volumes of approximately 18 and 32 
liters per min. Such curves may be used as a basis for design of 
cams for driving mechanical pumps. In order to obtain the data 
required for the construction of these curves, the minute volume 
must be predetermined. 
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Fic. 2. Maximum and minimum inspiratory flow-time curves 
obtained by tracing of original records. The minimum flow-time 
curve pertains to a subject with the lowest ventilation and the 
lowest respiratory rate in a group of 27 sedentary subjects. These 
respiratory curves may also be used for the design of test equip- 
ment for respiratory apparatus. 


under consideration, the value becomes 66.0 liters per 
min. B.T.P.S. The maximum figure actually observed 
in these flights was approximately 60 liters per min. 
B.T.P.S. 

In the design of equipment the average figures for 
ventilation have value for use in design of mechanical 
lungs for fatigue testing, for varying degrees of work, 
and also are utilized for mechanical mask freezing tests, 
etc. The ventilation must be further analyzed to deter- 
mine maximum, minimum, average, and instantaneous 
flow patterns (Figs. 1, 2). Data on these at present are 
taken from Silverman, ef al.’ (that data which involves 
no inspiratory resistance), since a complete study is 
available for 27 subjects. Studies of this nature have 
also been made by Hull, et a/.8 General deductions from 
these data have been transferred to data collected in 
aircraft. The minimum ventilation given is 6.1 liters 
per min. B.T.P.S.; the maximum is 90 liters per min. 
B.T.P.S. (Compare with minimum of 6.8 and maxi- 
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mum of 60 liters per min. in aircraft experiments.°) 
Average respiratory patterns in use for fatigue testing 
and mask freezing have been made as composites from 
the data of Lovelace, ef al.,5 and Silverman, et al.’ 
(Fig. 1). 

The effect of activity on ventilation has been well 
shown by Silverman, al.,’? and Wilson, et Earlier 
work has been reported by Dill’ and Penrod." The 
effect of cold would seem to be negligible from data of 
Lovelace, ef al.,5 although other data indicate an in- 
crease in ventilation when subjects were exposed to the 
cold. Flight conditions seem to effect an increase in 
the ventilation of fighter personnel with increasing alti- 
tude and of bomber personnel at altitudes above 33,000 
ft.2 Altitude introduces a factor above 33,000 ft. ™ 
The effect of clothing, tight-fitting parachutes, and 
position in the airplane on respiratory patterns has not 
been quantitatively defined. Dill has shown that pro- 
tective clothing may add somewhat to the oxygen 
requirement for performing a task, particularly if the 
clothing is ill fitting, and that exposure to cold (to the 
point of shivering) may increase the oxygen require- 
ment six times. The actual flight data of Carlson and 
Wulff!? seem to indicate that the factor of clothing and 
cold is not evident or is a constant in their overall 
averages. More recently Ferguson has further eval- 
uated this problem on a larger group of subjects and his 
results are shown in Table 1. Since it has been shown’ 
that the inspiratory pattern changes with the type of 
exercise, further work is indicated. 

Instantaneous rates of inspiratory gas flow are of 
interest from several points of view. The variation in 
flow during respiration must be used to duplicate respir- 
atory patterns. Maximum instantaneous inspiratory 
flow indicates the maximum flow to be met by the res- 
piratory equipment. The mean inspiratory flow has 
been used by Hart?* in mask leakage considerations. 

An analysis of the data of Silverman, ef al.,’ with the 
design of oxygen dispensing equipment specifically in 
mind discloses a number of relationships which may 
serve as a physiologic basis and form a framework which 
must be tested experimentally and upon which future 
data may be applied. These data furnish basic in- 
formation about four important respiratory phenomena: 
minute volume, respiratory rate, maximum instantane- 
ous flow, and portions of the respiratory cycle spent in 
inspiration and expiration. When the latter three are 
plotted against the minute volume certain regularities 
become apparent (Fig. 3). The relative time of in- 
spiration is less than half at low ventilations and in- 
creases to one-half with increasing ventilation. The 
maximum instantaneous flow may be determined from 
these data as given in Fig. 3 as: 


M.I.F. = 2.6V + 10 (4) 
The 90 per cent inclusive maximum instantaneous 
flow may be calculated from: 


M.L.F. (90) = 3.1V + 16 (5) 


DESIGN OF RESPIRATORY EQUIPMENT 


TABLE 1* 


Effect of Temperature, Anoxia, and Work on Pulmonary Ven- 
tilation, measured by inhaling air through a dry gas meter wear- 
ing R.C.A.F. oxygen mask. 

Number of subjects—66 (normal adult males). 

Work—Standing up (on floor) and sitting down (on bench) 12 
times per min. for 3 min. Volume measured for the third minute 
of exercise. 

Seated—Sitting at ease but not absolutely quietly. 

Cold——25° to —40° with wind 10 m.p.h. Subjects sat 15 
min. in cold before measurements taken, first seated then with 
exercise. 


Clothing—Flying clothing in cold. Indoor clothing at room 


temperature. 
G.L.—585 ft. 
Pulmonary Venti- 
Liters lation 
Altitude per Min. B.T.P.S. 
Temperature in Ft. Activity Range Average 
Room temperature G.L. Seated 4.3-20.7 10.2 
Room temperature G.L. Work 7.3-34.3 16.5 
Room temperature 10,000 Seated 5.1-19.9 11.0 
Room temperature 10,000 Work 9.5-33.3 19.5 
Room temperature 15,000 Seated 6.4-19.1 11.9 
Room temperature 15,000 Work 9.1-37.1 22.6 
Cold G.L. Seated 8.6-31.1 17.1 
Cold G.L. Work 7.7-54.6 28.3 
Cold 10,000 Seated 4.0-34.6 17.8 
Cold 10,000 Work 10.6-49.7 31.6 
Cold 15,000 Seated 8.1-24.8 17.1 
Cold 15,000 13.0-58.6 33.2 


Work 


Most obvious source of error is mask leakage which would tend 
to make results low. 

* Data from No. 1 Clinical Investigation Unit R.C.A.F. to be 
published. Furnished through the courtesy of J. K. W. Ferguson, 
M.D. 
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Fic. 3. Presentation of average respiratory data. This graph 
represents several factors in respiration plotted against the 
minute volume. The data used for all curves except 90 per cent 
inclusive maximum instantaneous flow (M.I.F.) curve are average 
data obtained from reference 7. 

The 90% inclusive M.I.F. curve was calculated from the aver- 
age minute volume data of reference 7. The standard error of the 
M.I.F. 90 per cent data, computed for the data series having the 
largest standard deviation, is +0.63 liters per min., indicating 
that out of any number of future groups of data, 90 per cent of 
the M.I.F. 90 per cent figures will fall within +2 liters per min. of 
the indicated values. 
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Fic. 4. Inspiratory flow-time pattern marked to obtain data 
for construction of curves on Fig. 5. An average inspiratory flow 
curve for a work load of 620 Kg. meter per min. obtained from 
reference 7. This curve is divided into seven parts by vertical 
lines which intersect the flow curve at 40, 80, and 120 liters per 
min. Relative flow is obtained from the ratio of fractional flow to 
maximum flow; relative volume from the ratio fractional time to 
total time. A plot of the latter two against relative flow is illus- 
trated in Fig. 5. ; 


The 90 per cent inclusive minute volume as taken from 
reference 7 may be determined by the equation: 
V (90) = 14V — 2 (6) 


The mean inspiratory flow is 0.9 the maximum in- 
stantaneous flow and is: 


mean I.F. = 2.3V + 9 (7) 
and 
mean I.F. (90) = 2.8V + 14 (8) 


Figs. 4-7 demonstrate that the actual inspiratory pat- 
terns can be generalized. The data from Figs. 20-24 
of Silverman, et al.,’ representing ventilations from 8 to 
40 liters per min. are shown plotted as fraction of tidal 
volume inspired versus flow as a fraction of the maxi- 
mum flow. The two figures supply information for 
determining representative curves at any minute 
volume. 

Expiratory rates of flow should be categorized simi- 
larly. Further data must be obtained. These data are 
important in design of expiratory valves. They cannot 
be transposed directly to intermittent or continuous 
pressure breathing equipment since the effects of the 
pressure alter the respiratory pattern.’ Silverman’s 
data (unpublished) indicate that the expiratory pattern 
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will differ from the inspiratory pattern. Maximum ex- 
piratory rates of flow are found to be 2.8 times the 
ventilation. 

The per cent of tidal volume drawn in at various rates 
of flow has proved significant in certain aspects of con- 
tinuous flow equipment, in the design of demand equip- 
ment, and in the consideration of mask leakage. Hart 
indicates that about 29 per cent of the inspired volume 
of air is inhaled at rates below 75 per cent of peak 
flow rates and about 60 per cent of the inspired volume 
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Fic. 5. Generalized relative volume—relative flow and rela- 
tive time—relative flow curves obtained in the manner described 
in the explanation for Fig. 4. From these curves average data 
may be extracted for construction of inspiratory flow curves and 
for construction of curves in Fig. 7. 
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Fic.6. Generalized relative volume-relative flow curve marked 
to illustrate extraction of data used to construct curves of Fig. 
7. The flow ratio is arbitrarily selected, i.e., 0.50, and a horizon- 
tal line intersecting the relative volume curve is drawn. The 
ratio of the area below this horizontal line and bounded by the 
relative volume curve (hatched and cross-hatched region) and 
the total area is plotted against relative flow to obtain data for the 
construction of curve 2, Fig. 7, labeled ‘‘Total Volume at Frac- 
tion.’’ Where the horizontal line, i.e., 0.50, intersects the rela- 
tive volume curve, perpendicular lines are dropped to the abscissa. 
The ratio of the sum of areas between these perpendicular lines 
and the relative volume curve (cross-hatched lines) and the total 
area beneath the relative volume curve is plotted against relative 
flow to obtain curve 1 of Fig. 7, labeled ‘‘Total Volume at Frac- 
tion or Less.” 
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is inhaled at rates below 95 per cent of the peak rate at 
rest and during exercise. This agrees with Fig. 7. 
The individual variation is great and the report states 
that during quiet breathing 1-14 per cent of the in- 
spired volume is inhafed at rates below 10 liters per 


min. 


RESPIRATORY RATE 


The general relationship of respiratory rate to venti- 
lation is shown in Fig. 3. At low ventilations, between 
5 and 20 liters per min., this curve may be sigmoid al- 
though a scatter diagram of the actual data indicates a 
straight line is more applicable. In this range the 
variation in respiratory rate is great. The data of 
Silverman’ indicate that where the average minute 
volume is 16.6 the range is from 8 to 24 breaths per min. 
The extent of the range at any given average ventila- 
tion decreases with increasing ventilation. 


HuMAN TOLERABLE RESISTANCE TO FLOW 


Tolerable inspiratory and expiratory resistance to 
flow has been indicated! and the effect of resistance on 
flow rates studied.” Negative pressures less than 1 in. 
water column are unnoticed at 50 liters per min. flow 
(see Fig. 8). At 100 liters per min. flow suctions as high 
as 1.5 in. water column are unnoticed. Maximum 
instantaneous rates of flow decrease and the relative 
time of inspiration increases with increased inspiratory 
resistance. Increasing inspiratory resistance decreases 
the minute volume and the maximum instantaneous 
rate of flow,and the relative timeof inspiration increases. 
At high work loads a resistance of 50 mm. (measured 
at 80 liters per min., 70°C. dry flow) causes a 7 per cent 
decrease in ventilation and approximately 15 per cent 
decrease in the maximum instantaneous flow. Demand 
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Fic. 7. Generalized fractional volume-relative flow curves. 
The curves in Fig. 7 present generalized inspiratory flow data ob- 
tained from curves in Fig. 6 as follows: (a) Curve 1 is a plot of 
the ratio of flow (instantaneous flow and maximum instantaneous 
flow) against the fraction of the volume inspired at a given in- 
stantaneous rate of flow. (b) Curve 2 is a plot of the ratio of the 


flow against total volume inspired at a given instantaneous rate of 
flow. These data are obtained from Fig. 6 in the manner described 
in that figure. 
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Fic. 8. Resistance to breathing. Data illustrating negative 
mask pressure in relation to maximum instantaneous rates of 
flow. The data were obtained by the instantaneous measure- 
ment and recording of pressure in a mask fitted with orifices of 
varying diameter. The open circles denote resistances which 
were not noticed by subjects. The closed circles denote resist- 
ances noticed by the subjects; the closed triangles denote re- 
sistance uncomfortable for short periods. Curve A represents the 
maximum allowable inspiratory resistance that is unnoticed by 
subjects; curve B represents the maximum allowable resistance 
that is noticed but not classified as uncomfortable by subjects. 
Data are from reference 1. 


and continuous flow equipment require low resistance or 
positive pressure and are well below tolerable limits 
generally. Expiratory resistance is necessary in pres- 
sure demand equipment. 


Gas EXCHANGE IN THE LUNGS 


For the design of closed respiratory systems it is 
necessary to consider the magnitudes of the gas ex- 
changes that take place in the lungs (Fig. 9). These 
changes are: (1) the removal of oxygen from inspired 
gas (oxygen consumption); (2) the addition of carbon 
dioxide to expired gas (CO, elimination); and the 
addition of water vapor to expired gas (H2,O vapor 
elimination). Oxygen consumption and carbon dioxide 
production vary among individuals and the respective 
curves are experimentally determined and represent 
the 90 per cent inclusive figures for oxygen consumption 
and carbon dioxide elimination. The data represent 
the volume of vapor in expired air saturated at 37°C. 
Seeley’ has published data on the rate of saturation of 
expired air and indicates that it approaches 100 per 
cent. Burch” has reported that expired air is not and 
cannot be saturated with water. The average relative 
humidity of expired air in a comfortable environment is 
reported to be 88 per cent at a temperature of 33°C. 


APPLICATION TO DESIGN 


The application of these data can be sampled by 
reference to the schematic types of continuous flow 
systems in Fig. 10. In a straight continuous flow 
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Fic. 9. Curves illustrating the relation of oxygen consump- 
tion, carbon-dioxide elimination, and water-vapor elimination to 
minute volume. The oxygen consumption and carbon-dioxide 
elimination curves are based on the 90 per cent inclusive figure 
calculated from average data obtained from Silverman, un- 
published data. The reliability of these curves was checked by 
plotting data from Schneider, Physiology of Muscular Activity. 
Open squares represent carbon-dioxide elimination and closed 
squares represent oxygen consumption. Notice that these values 
taken from Schneider fall on or below the 90 per cent inclusive 
curves. The standard error of the oxygen consumption and car- 
bon-dioxide elimination 90 per cent inclusive values, computed 
fcr the data series having the largest standard deviation, is 
+0.0047 liter per min., indicating that out of any number of 
future groups of data, 90 per cent of the oxygen consumption and 
carbon-dioxide 90 per cent inclusive values will fall within 
+0.014 liter per min. of the indicated values. The curve for 
water vapor elimination is calculated, assuming the expired air is 
saturated with water vapor at 37°C. See text for further discus- 
sion. 


system (Fig. 10-1), such as the Japanese use, the con- 
tinuous supply of oxygen must be carefully related to 
the volumes inspired under given instantaneous flows. 
A use of the data for straight continuous flow may be 
given. Since only that part of the flow is used which is 
flowing during inspiration, the flow required will be 
equal to the flow at which a fraction of the inspired 
volume equals the fraction of oxygen required. Thus, 


F = M.I.F. X fraction of M.I.F. under which frac- 
tion of tidal volume taken in is equal to the 
fraction of added oxygen required. (See 
Fig. 7.) 

since: 
M.1.F. = 2.6V + 10 (9) 


Fraction added oxygen required = 1.27 X 
(fO2 — 0.21) (10) 


(converting oxygen required to fraction of oxygen 
added). Fraction of M.I.F. under which fraction of 
tidal volume is taken in is equal to the fraction of added 
oxygen required taken from Fig. 4, Curve 2 and be- 
tween the fractions 0.1 and 0.9 is: 


(8/9) [1.27(fO2 — 0.21)] (11) 


The flow required for a continuous flow system without 
economizer systems is: 


F = (2.6V + 10){(8/9)[1.27(fO2 — 0.21)]} x 
{[(P. — 47)/P](T/T,)} (12) 


where the fraction {[(P, — 47)/P](T/T,)} corrects 
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the values to 760 mm. Hg and 0°C. This equation, be. 
cause of the inadequacies of the fraction 8/9, is valid 
when fO: is greater than 0.10 and less than 0.90. 

These equations may be simplified by the substitu- 
tion of the actual fraction of oxygen desired. 

Flows calculated from these data are in fair agree- 
ment with Swann.'® When a reservoir or rebreather 
reservoir (Fig. 10-2 and -3) are used the relationship 
of the factors change. These systems are more complex 
to evaluate than the demand system. Making the 
assumption that the reservoir will save half the oxygen 
the above equation can be utilized. Values here are 
analogous to those reported by Swann. Further 
refinements of the formula to evaluate the effect of a 
rebreather reservoir can be made from the graphic data. 
Limitations of the rebreather economizer due to its 
volume can be further expressed by: 


(V/n)/V, (13) 
and the equation becomes: 
2F = (2.6V + 10){ (8/9) [1.27(fO2 — 0.21)]} x 
(Pa — (14) 
This equation is in error by the amount of flow that oc- 


curs during the inspiratory cycle and can be further cor- 
rected with the use of: 


{[(V/n)/V,] — (8/9)[1.27(fO2 — 0.21)]}/2 (15) 


An economizer without a rebreather will require the 
modification of the basic formula whenever [(V/n)/V,] 
is greater than one. With a 700-cc. economizer, 
[(V/n)/V,] becomes 1.0 at 15 liters per min. ventilation 
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Fic. 10. Schematic sketches of oxygen supply systems are 
illustrated. Sketch 1 represents a constant flow system consist- 
ing only of regulator and mask. Sketch 2 represents a constant 
flow system consisting of a regulator, mask, and a reservoir 
volume. Sketch 3 represents a constant flow system consisting 
of a regulator, mask, and a rebreather bag. Sketch 4 represents 


a demand system consisting of a regulator and mask. Refer to 


text for further detail. 
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which may lead to the statement that little difference 
exists in the two systems.” Since these formulations 
incorporate the fraction of oxygen required they may be 
utilized to calculate flows necessary for oxygen therapy 
and such calculations show good agreement with 
Boothby*! and unpublished data from the Aero Medical 
Laboratory. 

In the demand system, requirements may be taken 
from the Figs. 3, 5, and 7 when the value of pulmonary 
ventilation is set. Army Air Forces’ equipment is de- 
signed for an average inactive ventilation of 14 liters per 
min. B.T.P.S.; an average active ventilation of 24 
liters per min. B.T.P.S.; and a maximal ventilation of 
60 liters per min. B.T.P.S. Using Eq. (6) the 90 per cent 
inclusive ventilations are 18 and 32 liters per min. 
B.T.P.S., respectively. Using Eq. (4) the 90 per 
cent inclusive maximum instantaneous flows are 45, 
93, and 166 liters per min., respectively. Oxygen 
regulators must be designed to deliver these in- 
cdusive flows within the minimum suction require- 
ments. 

In the testing of respiratory equipment, average data 
should be employed and corrected to the 90 per cent 
inclusive value. At the present time A.A.F. regulators 
are tested with a mechanical lung that duplicates the 
breathing pattern for 18 and 32 liters per min. ventila- 
tion, illustrated in Fig. 1. These flow-time curves were 
reconstructed from the curves in Figs. 3 and 5. In addi- 
tion to testing with average +1.3 standard deviation 
data, extremes of the range encountered should also be 
used. Thus, in the inactive ventilation group’ there is 
illustrated the flow-time pattern of the lowest ventila- 
tion (Fig. 2), that should be duplicated on the mechani- 
cal lung. A ventilation pattern for 70 liters per min. 
(Fig. 2) is used to test performance at maximum respira- 
tory volumes. 

The variables of respiratory phenomena discussed 
above have not been clearly defined for pressure breath- 
ing equipment.’ Interrelated physiologic factors as 
well as characteristics of equipment complicate the 
problem but do not obscure the solution. A similar 
approach can and will be made to reduce these factors 
to generalities for equipment design. 


CONCLUSION 


The synthesis of varied and apparently independent 
data into a unified and generalized relation is desirable 
when large and diverse groups of data exist and this 
data is to be transferred from the field of physiologic 
investigation to serve as a basic requirement for an 
engineering development. In addition, such a syn- 


thesis has the advantage of establishing a measure for 
evaluation of such data, readily directs research into 
channels where sufficient data are lacking, and permits 
tapid and sufficiently accurate (i.e., 10 per cent) ex- 
traction of average data for practical application. The 
presentation of the general relations of respiratory data 


in the foregoing report embodies the characteristics of 
such a synthesis. 

The construction of such a generalized framework is 
based upon experimental data. The value of such a 
framework depends upon the accuracy and validity of 
the data used. As an example, the average curve for 
respiratory rate in Fig. 3 embodies a large range among 
individual determinations which is not apparent upon 
inspection of the curve but is apparent when compared 
with other experimental data. Similarly, the curve for 
water vapor elimination in Fig. 9 is not based on experi- 
mental data but is calculated on the basis of certain 
assumptions (refer to text). Comparison of this curve 
with experimental data, as these become available, 
shows a discrepancy, since the values in Fig. 9 are 
probably too high. 

In spite of the faults of the framework presented in 
this report, it is believed to be a worth-while approach 
to a synthesis of respiratory data which should not be 
considered a definitive attempt but rather a directional 
attempt to be modified and corrected by subsequent 
experimental data. 
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Letter to 


Dear Sir: 

As is pointed out in a Letter to the Editor in the JoURNAL OF 
THE AERONAUTICAL SCIENCES, Vol. 14, No. 1, p. 23, January, 
1947, by Messrs. Krzywoblocki and McCloy, a considerable 
amount of confusion exists with regard to the concept of pro- 
pulsive efficiency. This quantity, 7), is defined as the ratio of 
“useful work” or product of thrust and air speed, to “‘energy 
input” to the gases ejected from the propulsor. For the con- 
ventional propeller, it reduces to 


np = 2/[1 + (V,/V)] (1) 
where V, is the air speed in the slipstream and V is the aircraft 
speed. For a rocket with exhaust gas ejection velocity C, np 
reduces to 

np = 2V/C (2) 
a quantity that varies from 0 to ~ as V increases without 
limit. 

The confusion that exists with regard to y, is blamed by the 
authors of the previous letter, on the fact that this efficiency can 
have a value greater than 1. It would seem that there is more 
wrong with yn, than that. If, for example, the expression for V 
as reached by a rocket in empty space, is put into Eq. (2) above, 
the result is disturbing. For such a rocket, 


V = log.(mo/m) (3) 
where c¢ is the exhaust velocity of the combustion products, mo 
and m are masses of the craft when V is 0 and V, respectively. 
It follows from Eq. (3) that f 

np = 2 loge (m/m) (4) 


The “‘efficiency”’ is here independent of the fuel and even of the 
propulsive mechanism. The change in definition of , proposed 


the Editor 


by the authors mentioned does not improve its properties in this 
respect. 

It is suggested that the whole basis of definition of », should 
be re-examined. The idea of an efficiency 7 is useful in certain 
fields of engineering. In particular, it characterizes a machine 
that does useful work. When 7 is 100 per cent, there are no 
friction losses; » cannot be greater than 100 per cent without 
violating the energy conservation law. 

If a quantity 7, is applied to the flight process, a great difficulty 
is at once apparent: there is no standard ‘‘useful work’’ required 
to accomplish flight that cannot in principle be reduced toward 
0 as a limit. If the thrust power necessary to be expended in 
steady flight is defined as ‘‘useful work,’”’ any value of n, > 100 
per cent does not indicate a failure of the conservation of energy; 
if anything, it simply indicates an error in definition of the 
“useful work”’ required to fly. 

Furthermore, thrust power devoted to acceleration rather 
than overcoming drag, is not ‘‘expended”’ in the ordinary sense 
of the term as used in defining »,. Such power is really stored in 
the form of kinetic energy of the craft. It can be applied to the 
overcoming of drag in flight at a later time. Therefore, the 
use of 7, to characterize accelerated flight in any sense at all— 
and this is especially true for the case of rockets—is bound to 
lead to confusion and is probably meaningless. 

It is suggested that specific impulse in the form c, or J = ¢/g, 
where g is 32.2 ft. per sec.?, is a proper performance criterion for 
rockets. The impulse per unit of fuel mass expended is proposed 
as preferable to , also in characterizing propeller-driven aircraft. 

Davip T. WILLIAMS 
Associate Professor of Aeronautical 
Engineering 
University of Michigan 
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Numerical-Graphical Methods of 
Characteristics for Plane Potential 


Shock-Free Flow Problems 


LESTER L. CRONVICHt 
Applied Physics Laboratory, The Johns Hopkins Unwersity 


SUMMARY 


Two numerical-graphical methods of solution of supersonic 
flow problems are developed from the basic theory of charac- 
teristics. Both are restricted to plane potential shock-free 
flow. The ‘‘lattice-point’’ method, recommended for flows with 
given boundaries, is applied to flow in a hyperbolic channel. 
The “‘field’’ method, suitable for problems in which flow bound- 
aries are to be determined, as in the case of supersonic nozzle 
design, is used for a problem of parallel routing of an initially 
expanded flow. Complete rules of procedure are presented for 
each method. 


INTRODUCTION 


N 1929 PRANDTL AND BUSEMANN! presented a 
| graphical scheme for determining two-dimensional 
supersonic flows, which has come to be known as the 
Prandtl-Busemann ‘“‘method of characteristics.’’ The 
essential features of the method are explained in Buse- 
mann’s handbook article? and more recently in the 
book by Sauer.* A detailed exposition of the funda- 
mental ideas and procedures is contained in a report 
by Heybey.* 

The mathematical basis for the method may be found 
in the book Methoden der Mathematischen Physik, by 
Courant and Hilbert’ or in a series of reports by Toll- 
mien.’ A brief discussion of the theory necessary for 
the methods of solution of this paper is given in the 
following section. 


MATHEMATICAL THEORY 
Consider a differential equation of the form 
Ado: + Bory + Chy +D=0 (1) 


where A, B, C, and D are functions of x, y, ¢, ¢,, and 
¢, and satisfy the relation 


— 4AC >0 


The function ¢(x, y) may be considered to define a 
surface in xyz-space. Certain space curves lying on 
this surface have projections in the xy-plane which 
satisfy the differential equation 


Ady? — Bdydx + Cdx? = 0 (2) 
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These projections are usually referred to as ‘‘charac- 
teristics.’’ The properties of these curves® make 
possible the development of the approximate methods 
discussed here. From Eq. (2) it is found that there 
are two families of characteristics in the wxy-plane 
given by 


dy/dx = (B — VB? — 4AC)/2A 
dy/dx = (B + VB? — 4AC)/2A 


One of the properties of the characteristics is that 
the differential equation, Eq. (1), can be satisfied along 
such a curve without having to depart from it obliquely. 
In other words, if initial values of ¢, ¢,, and @, are 
given on a curve C which is intersected by characteris- 
tics, one may proceed along characteristics and de- 
termine the unknown values of ¢, ¢,, and ¢, on the 
curves by a step-by-step process in the xy-plane. Since 
the region of definition of @ in the xy-plane can be 
covered by a mesh of characteristics of any desired 
fineness, it is evident that theoretically a solution 
may be obtained with any degree of approximation 
desired. In practice, however, the degree of approxi- 
mation is limited by the amount of numerical or 
graphical work involved in obtaining such a solution. 
The idea of solving the differential equation by pro- 
ceeding along characteristics is the basic idea in the 
development of the ‘‘lattice-point’’ method of solution. 

Another property of characteristics is that they are 
curves across which the second derivatives of @ may 
experience discontinuities. Thus, solutions of differen- 
tial equation, Eq. (1), in two regions separated by a 
characteristic may be patched together along this 
characteristic. This is the basic idea in the develop- 
ment of the “‘field’’ method of solution of this paper 
and also that of Prandtl and Busemann. 


(3) 


PLANE, SHOCK-FREE, POTENTIAL FLOW 


For plane, shock-free, potential flow the differential 
equation of motion is 


2uv v? 


where 


= velocity potential 
u = ¢, = x-component of velocity 


7 


° 

= 
| 


v = gy = y-component of velocity 

a = V[(k — 1)/2](qm? — = local velocity of 
sound 

= Vu? + v? = velocity 

dm = maximum velocity, corresponding to a = 0 

k = ratio of specific heats 


Comparing this equation with Eq. (1) it is found that 
D = Oand 


B* — 4AC = 4[(qg?/a?) — 1] > 0 


for g > a, that is, for supersonic velocity. 
The two families of characteristics are 


dy wtaVq—a? 


dx u? — a? (5) 
5 

dy w—aVq—a? 

dx u? — a? 


Now let @ be the inclination of the velocity vector 
with the x-axis, so that 


v =qsin#@ (6) 
By differentiating, 


du = dqcos @ — qdé sin 0 
dv dq sin 6 + qd@ cos 6 


Now at any point P let the x-axis coincide with the 
stream direction so that Eqs. (6) and (7) become 


u=q, v=0 (6a) 
du = dq, dv = qdé (7a) 


(7) 


Using these transformations the differential equation, 
Eq. (4), may be written 


[1 — (q*/a*)]dg(dy/dx) + gd = 0 (8) 


Furthermore, the equations of the characteristics be 
come 


(9) 


where a@ is the well-known Mach angle usually defined 
by the relation 


sina = a/qg = 1/M (10) 


M being the Mach Number. 

Thus one family of characteristics forms the Mach 
angle a with the tangents to the streamlines and is 
referred to as “‘left-running’’ characteristics or Mach 
lines, since any one of them is on the left as one faces 
in the downstream direction; the other family forms 
the angle— a with the tangents to the streamlines and 
is referred to as “right-running”’ characteristics or 
Mach lines. 
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Fic. 1. Pressure ratio p/p), Mach Number M, and Mach angle 
a vs. v for plane potential flow. 


The differential equation, Eq. (8), can now he 
written with the aid of Eqs. (9) in the form 
dq — (a/Vq — a*)qdo = 0 

dq + (a/V — a*)qd8 = 0 
along the left-running and right-running Mach lines, 


respectively. These equations can be integrated to 
yield 


(11) 


6 = v(q/a) + Ky, 
where 
are tan (13) 
a 


K,, Kr = constants of integration. One method for 
carrying out this integration is indicated on p. 56 of 
reference 3. The results presented there can be trans- 
formed to Eqs. (12) above. 

By application of Bernoulli’s equation it is possible 
to express the pressure ratio p/po, where po is the 
pressure corresponding to zero velocity, as a function 


of g/a. Thus 
2\ —k/(k — 1) 
42512) (14) 


The Mach angle a may likewise be expressed as a 
function of g/a by the equation 


a = arc sin (a/q) (10a) 


In view of Eq. (13), therefore, Mach Number M = 
q/a, Mach angle a, and pressure ratio p/p) may be 
regarded as functions of v. For convenience in ap- 
plication these quantities are plotted against ». A 
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Fic. 2. Construction of characteristics for flow in a hyperbolic 
channel, M = 1.6. 


sample of such curves for air (k = 1.405) is given in 
Fig. 1. 


“LATTICE-PoINT’”? METHOD OF SOLUTION 


In the preceding section the characteristics were 
identified with the familiar Mach lines in the flow plane. 
These characteristics form a mesh covering the flow 
region, and their points of intersection give a system 
of “lattice points.” In the “‘lattice-point’’ method of 
solution one proceeds along the characteristics from 
lattice-point to lattice-point (the characteristics being 
made up of straight-line segments for this approxi- 
mation). All calculations along left-running charac- 
teristics are governed by the equations 


dy/dx = tana \ 
6 = v(g/a) + 


and along the right-running characteristics by the 
equations 


(15) 


dy/dx tan (— a) (16) 


= —v(q/a) + Kr 
The method is to be used for flows with given boundaries 


10 
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Fic. 3. Angle of inclination of channel wall of Fig. 2 in degrees. 


for which the Prandtl-Busemann method and the 
“field’’ method of this paper are not always suitable. 

The general procedure of this method can be given 
by the following set of rules which are applied to flow 
in a hyperbolic channel with initially uniform flow at 
a Mach Number of 1.6 (see Figs. 2 and 3 and Table 1). 

(1) Construct a scale drawing of the flow boundaries 
(as in Fig. 2). 

(2) Plot 6,, angle of inclination of the wall to the 
initial flow direction. This is illustrated in Fig. 3. 

(3) Choose a suitable set of points on the initial 
curve (points 0, 1, and 2 in Fig. 2), and from the given 
initial conditions determine v for each of the points 
from Fig. 1. (These are recorded in Table 1.) Points 
of intersection of the initial curve and the boundaries 


TABLE 1 
Calculated Data for Hyperbolic Channel of Fig. 2 
R L L R 

Pt @° (@+»)° (@—»)° a (@+a)° (—a)° M b/po x 
0 0 14.9 14.9 —14.9 38.6 38.6 —38.6 1.6 0.234 0 

1 0 14.9 14.9 —14.9 38.6 38.6 —38.6 1.6 0.234 0 

2 0 14.9 14.9 —14.9 38.6 38.6 —38.6 1.6 0.234 0 

3 0 14.9 14.9 —14.9 38.6 38.6 —38.6 1.6 0.234 0.63 
4 0 14.9 14.9 —14.9 38.6 38.6 —38.6 1.6 0.234 0.63 
5 —2.6 17.5 14.9 —20.1 36.2 33.6 —38.8 1.69 0.205 1.31 
6 0 14.9 14.9 —14.9 38.6 38.6 —38.6 1.6 0.234 1.26 
7 2.6 17.5 20.1 —14.9 36.2 38.8 —33.6 1.69 0.205 1.31 
8 —2.6 17.5 14.9 —20.i 36.2 33.6 —38.8 1.69 0.205 1.98 
9 2.6 17.5 20.1 —14.9 36.2 38.6 —33.6 1.69 0.205 1.98 
10 —5.0 19.9 14.9 —24.9 34.3 29.3 —39.3 1.775 0.181 2.67 
ll 0 20.1 20.1 —20.1 34.1 34.1 —34.1 1.78 0.179 2.87 
12 5.0 19.9 24.9 —14.9 34.3 39.3 —29.3 1.775 0.181 2.67 
13 —2.4 22.5 20.1 —24.9 32.4 30.0 —34.8 1.825 0.157 3.67 
14 2.4 22.5 24.9 —20.1 32.4 34.8 —30.0 1.825 0.157 3.67 
15 -—7.5 27.6 20.1 —35.1 29.2 21.7 —36.7 2.05 0.1185 4.85 
16 0 24.9 24.9 —24.9 30.8 30.8 —30.8 1.95 0.138 4.65 
17 7.5 27.6 35.1 —20.1 29.2 36.7 —21.7 2.05 0.1185 4.85 
18 —5.1 30.0 24.9 —35.1 27.9 22.8 —33.0 2.135 0.103 6.10 
19 5.1 30.0 35.1 —24.9 27.9 33.0 —22.8 2.135 0.103 6.10 
20 —9.0 33.9 24.9 —42.9 25.9 16.9 —34.9 2.29 0.081 7.45 
21 0 35.1 35.1 —35.1 25.3 25.3 —25.3 2.34 0.075 8.16 
22 9.0 33.9 42.9 —24.9 25.9 34.9 —16.9 2.29 0.081 7.45 
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should be included. It is recommended that one 
make use of the symmetry of the figure as a check on 
the accuracy of construction. If a discontinuous 
change of 0, occurs at the intersection of the initial 
line and the boundary, several coincident points should 
be taken, each of which accounts for a fraction of the 
total jump. The coincident points lie on the same 
left-running characteristic if they are on the upper 
boundary (flow direction being from left to right) and 
on the same right-running characteristic if they are 
on the lower boundary. 

(4) Draw both left-running and right-running 
characteristics (straight-line segments) at each point 
chosen on the initial line provided these lines extend 
into the flow region. The inclinations of these line 
segments to the initial flow direction are given by the 
columns (@ + a) and (@ — a) of Table 1. Their 
intersections locate new lattice points on the charac- 
teristics. 

(5) Values of 6 and » for each new point are ob- 
tained by simultaneous solution of Eqs. (12). For 
example, for point 3 


0s — Vy = 
since points 0 and 3 lie on a left-running characteristic 
and 

=A +n 


since points 1 and 3 lie on a right-running characteristic. 
Thus we have 


+ v1) + (0 — w)]/2 
+ v1) —. (0 — vo)]/2 


(6) From Fig. 1 read off the value of a correspond- 
ing to vy at the newly found point. Calculate angles 
6 + a@ and @ — a@ and draw the characteristics as 
indicated in Step 4. 

(7) The procedure is modified in case one of the 
new points lies on the boundary (such as points 5, 7, 
10, etc., of Fig. 2). In such a case, the value of 6 
is read from the curve constructed in Step 2 (that is, 
Fig. 3 for this-illustration). And v is then obtained 
from whichever of Eqs. (12) applies. 

(8) Values of Mach Number and pressure ratio 
p/po at any point can be read from Fig. 1 after the value 
of v has been determined for that point. 


V3 


“FYELD’ METHOD OF SOLUTION 


As pointed out in the preceding section, the charac- 
teristics form a mesh covering the flow region, so 
that the region is broken up into small fields. In the 
field method of solution the characteristics are used 
for patching together the flows in two adjacent fields, 
throughout each of which the flow variables are as- 
sumed constant. Thus the crossing of a characteristic 
must account for a discontinuity or jump in certain 
flow variables. It is found convenient in each problem 


3 7 
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Fic. 4. Construction of characteristics and boundary for paral- 
lel routing of an initially expanded flow, M = 1.4. 


to assign a standard jump to the quantities 6 + » 
and @—v. The jump used for both quantities is 
numerically the same but opposite in sign. In ad- 
vancing from one field to the next a characteristic is 
crossed, thereby resulting in a jump in @ + » if a right- 
running characteristic is crossed or in 6 — » if a 
left-running characteristic is crossed. However, it is 
considered that the crossing of a right-running charac- 
teristic takes place by proceeding along an intersecting 
left-running characteristic (usually not drawn) along 
which the first of Eqs. (12) must be satisfied. Similarly, 
the crossing of a left-running characteristic takes place 
by proceeding along an intersecting right-running 
characteristic along which the second of Eqs. (12) must 
be satisfied. 

The slope of the boundary between two fields is 
determined from the values of @ and y in each of them 
in such a way as to satisfy approximately Eqs. (9). 
For example, if the boundary separating fields 1 and 3 
is a right-running characteristic, it is drawn so as to 
make an angle 6,. — a with the x-axis where 


= (A; + 63) /2 
and a@ is the Mach angle corresponding to vg, where 
Var = (1 + v3) /2 


Since this method, like that of Prandtl and Buse- 
mann, is suitable for problems in which flow boundaries 
are to be determined (such as free jets, nozzle design, 
etc.), the rules of procedure are laid down with this 
type of problem in mind. An application of the 
method is given in Fig. 4 and in Tables 2a and 2b. A 
steady flow at Mach Number 1.4, corresponding to 
y = 9°, enters a channel that initially turns the flow 
through 8°, making » = 17°, which corresponds to a 
Mach Number of 1.675. A channel must be so con- 
structed as to again direct the flow parallel to its 
original direction. 

The general procedure is as follows: 

(1) At the point of intersection of the initial line 
and the boundary (unknown beyond this point), 
determine v before and after the flow has been turned, 
from Fig. 1. The difference gives the total bend of 
the flow direction, vy — v, which is equivalent to 
6» — 4 of the upper boundary (the flow being from 
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TABLE 2a 
Field Data for Fig. 4 
‘Field 0° @—»)° M b/ bo 
0 0 9 9 -9 1.4 0.313 
1 4 13 17 -9 1.535 0.257 
3 8 17 25 -9 1.675 0.210 
2 —4 13 9 —17 1.535 0.257 
4 —-8 17 9 —25 1.675 0.210 
5 0 17 17 —-17 1.675 0.210 
6 -4 21 17 —25 1.815 0.170 
7 a 21 25 —17 1.815 0.170 
8 0 25 25 —25 1.955 0.137 
TABLE 2b 
Boundary Data for Fig. 4 
Boundary Var” a® (Oav + a) — a) 
0-1 2 11 42.9 —40.9 
0-2 —2 11 42.9 40.9 
1-3 6 15 38.5 —32.5 
244 —6 15 38.5 32.5 
1-5 2 15 38.5 40.5 
2-5 —2 15 38.5 —40.5 
4-6 —6 19 35.0 —41.0 
3-7 6 19 35.0 41.0 
5-6 —2 19 35.0 33.0 
5-7 2 19 35.0 —33.0 
6-8 —2 23 32.1 —34.1 
78 2 23 32.1 34.1 


left to right). The subscript 0 refers to initial values 
and F to values after the flow has been turned. In 
the illustration this angle is 8°. 

(2) Let m be the number of right-running charac- 
teristics to be crossed in changing from % + w to 
ip + vp along the upper boundary. The choice will 
depend on the accuracy desired. In Fig. 4, m = 2. 
The standard jumps are given by 


Ar = + vr — (00 + vo) |/n 
and 
A, = — Ar 


for the right-running and left-running characteristics, 
respectively. In the illustration the corresponding 
jumps are 8° and —8°. For a free jet the reflection 
of a characteristic in a boundary reverses the sign of 
the jump. Thus a reflected right-running character- 
istic with an original jump Ag becomes a left-running 
characteristic with a jump —A, = Ag. _ In other 
words, for a free jet each characteristic retains its 
original jump throughout the flow even though its 
direction is changed as a result of reflections. On the 
other hand, no change is encountered when a charac- 
teristic is reflected in a straight boundary. The 
reflected characteristic may be omitted provided the 
boundary is bent sufficiently to compensate for the 
omitted jump. This property is essential for the 
solution of the problem illustrated by Fig. 4. 

(3) Set up a table of values of 6, v, 6 + v, and 
4 — y for the various fields, starting with the initial 
values of 4 and v. To proceed to an adjoining field 


determine the type of characteristic that must be 


crossed. If it is a right-running characteristic, then 
in the new field 


6+ v = + w+ Ar 
and 
=  — w 
If it is a left-running characteristic, then 
6+ v = + w 
and 
— w+ Az 


From these pairs of formulas 6 and v can be found. 
The same procedure is followed in progressing from 
any known field to an adjacent unknown field. (See 
Table 2a.) 

(4) Set up a table of values of O5, var, a, Pas + a, 
and 6,, — a@ for the field boundaries as in Table 2b, 
using the procedure outlined earlier in this section. 

(5) All graphical construction can be given in this 
step since it amounts merely to the drawing of the 
boundaries between fields which make angles @,, + a 
and 6,4, — a with the original flow direction. Steps 
3, 4, and 5 are carried out simultaneously. 


COMMENTS ON THE METHODS 


The *‘Lattice-Point’’ Method 


Since the method is approximate, the effect of a 
change of slope of the boundary is carried along the 
flow only by the Mach lines being reflected from the 
walls. Consequently, the more Mach lines used, the 
greater is the possibility of showing the effects of the 
true contour. If the slope of the boundary is dis- 
continuous, a judicious selection of points along the 
initial line is needed so that the reflections of the Mach 
lines will occur as near to the points of discontinuity 
as possible. Brief preliminary estimates should be 
made to insure this behavior in the final solution. 

Theoretically the method is adaptable to problems 
in which flow boundaries are to be determined, as 
well as those with given boundaries. However, it is 
found in the case of the problem of a free jet expanding 
into air at lower pressure that the method of approxi- 
mation is such that the theoretical periodicity of such 
solutions is not obtained except in the limit—i.e., 
when infinitely many characteristics are used. This 
periodicity is obtained, however, in the solution by the . 
field method of this paper. 

One of the errors of approximation appears in the 
procedure of locating new points by drawing charac- 
teristics from known points with slopes given by data 
at these points rather than at some intermediate 
points. An iterative procedure may be employed if 
it is felt that the added computations are not pro- 
hibitive. Some improvement may be obtained in this 


. 
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way for flow regions in which the Mach Number is 
near unity. 


The ‘‘Field’’ Method 


The field method of this paper, although theoretically 
equivalent to the Prandtl-Busemann method, has an 
advantage over the latter in that jumps of amy magni- 
tude (fractional or integral) may be used. The Prandtl- 
Busemann method uses integral jumps necessitating 
the approximation of the problem itself as well as the 
solution. 

A disadvantage of both field methods is that they 
have not yet been suitably adapted to problems with 
arbitrarily given boundaries, such as in Fig. 2. The 
reason for this is basic. Standard jumps in 6 + »v 
and @ — »y are selected and must not be changed 
throughout the problem. Having chosen the jumps 
the directions of the Mach lines are then governed 
by the initial Mach Number. Thus one has no control 
over the locations of points of intersection of one 
boundary with Mach lines reflected from the opposite 
boundary. Consequently, it is highly improbable that 
any fields adjacent to the boundary would have a 
direction @ equivalent to the average direction of 
the section of the boundary forming one side of the 
field. 


The curves of Fig. 1 may be extended to higher Mach 
Numbers by use of Eqs. (13), (10a), and (14). Like. 
wise, they may be extended in the other direction to 
q/a = 1. For higher Mach Numbers greater care 
must be exercised in the graphical construction be- 
cause of the low values of 0+ aand@— a. Further. 
more, a different set of curves must be prepared for 
each value of k. The use of these curves rather than a 
table of pressure numbers is one of the principal ad- 
vantages of the method. 
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A Method of Airplane Performance 
Calculation Applicable to Any Polar 


R. E. WENDT* 


Western Air Lines 


SUMMARY 


A method is presented for the determination of airplane per- 
formance by use of formulas and charts. The conventional 
method of approximating a polar by expressing the effective in- 
duced drag as C,?/rA.R., is often inadequate because of the fact 
that many polars are not symmetrical about zero lift. A mathe- 
matical expression for the airplane polar is determined here by 
expressing the effective induced drag as a parabola symmetrical 
about the lift coefficient at which the minimum drag coefficient 
occurs. This coefficient is defined as C_,,;,, The minimum drag 
coefficient is defined as Cp,,;,, Hence, the airplane drag coeffi- 
cient is 


Co = Comin. + — 


Performance is determined by use of the fundamental loading 
parameters in common use, dp, As, and Az. 

The velocity at maximum lift drag ratio is used throughout as 
the basic velocity. All other velocities are expressed in terms of 
the basic velocity by use of the nondimensional velocity ratio, Ry. 

Three basic parameters are introduced: basic performance 
parameter, 

coefficient of asymmetry, 
Ca = Cimin./V A-R.eCoy 
performance parameter ratio, 
R,? 1 CaRy 


= — = 


Cp, is the equivalent parasite drag coefficient. The basic per- 


formance equation based on the usual assumption that lift = 
weight is 


Ai(dh/dt) = 1 — y’Ry 


When the effective induced drag term of the drag equation is 
expanded, a linear function of Cz appears. This term is ac- 
counted for by the coefficient of asymmetry. This coefficient has 
a constant value for any given polar. Rate of climb may be 
determined at any speed by use of the performance parameter Ry, 
which is a function of R,. Level flight velocity is determined in 
terms of R, for a given value of y and Ca. 

Because of the ability to approximate any polar with negligible 
error and the rigorous analysis used thereafter, the results 
obtained have a high degree of accuracy and are easily determined. 


NOMENCLATURE 


constant to convert horsepower to desired physical 
units 

A.R. = aspect ratio 

A.R., = effective aspect ratio=(A.R.), 

b = span 
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b.hp. = brake horsepower 

e = airplane efficiency factor 

f = equivalent parasite area = Cp;S 
dh/dt = rate of climb 

dh/dt; = indicated rate of climb 

L/D = lift-drag ratio 


= thrust 

Ts = time required to climb 

V = velocity 

VY = basic velocity or speed at maximum L/D 

W = airplane weight 

Ry = ratio of velocity to basic velocity 

Ry = performance parameter ratio = (R,*/2) + (1/2R,) — 
(CaRv/ Vx) 

Re = weight ratio = ratio of any weight to originally 
assumed weight 

Ap = parasite loading = 2W/pf loading 

rs = span loading = 2W/xpob*e parameters, 

Ne = thrust power loading = W/Ab.hp.m| physical 

units 

iy = parasite loading W/f loading param- 

a = span loading W/b?e eters, engineer- 

l, = thrust power loading W/b.hp.no! ing units 

Cdmin. = Minimum drag coefficient 

Ca = coefficient of asymmetry = Cy,,;,./ VAR.Co, 
Coa 

Cpa = asymmetrical drag coefficient = Cz,,;,./mA.R. 

Cpimin. = Minimum induced drag coefficient 

Crmin? = = Cra Cimin. 

Cp, = equivalent parasite drag coefficient = Cdmin, + 
Cdimin. 

Cimin. = lift coefficient at minimum drag coefficient 

y’ = basic performance parameter = ,d,°/4/d,)/4 = 

y = basic performance parameter (engineering units) = 
= 

vi = indicated performance parameter = po/p 

VR = basic climb performance parameter = YRy 

= basic velocity parameter = 

T = basic velocity parameter (engineering units) = (J,/,) 

A’ = fundamental performance parameter =),A;/(Ap/d.)* /* 

A = fundamental performance parameter (engineering 


units) = 
= propeller efficiency 


7 

0 = angle between flight path and horizon 

p = mass density of air 

Po = mass density of air at sea level 

= p/po 

\s» = Span-parasite loading = [po! 
(physical units) 

= span-parasite loading = W,/ (6Ve/f'/*) (engineering 
units) 


P; = b.hp.n Vo/(W/W,)*/? 

Veiw = Ve/(W/W.)*/ 

Ciw (dh/dt) Vo /(W/W,)* 

W, = standard weight or an assumed weight 
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Ryiw = Piw ratio = Piw/dep?/? 
AP = power in excess of power required for level flight 


INTRODUCTION 


N ANALYTICAL METHOD is presented here for the 
determination of airplane performance from test 
data or an estimated polar. The method is particularly 
applicable to polars determined from wind tunnel or 
flight test. Analytical methods of performance de- 
termination in common use express the effective induced 
drag simply as a function of C,?. It has been found 
that many test polars are not symmetrical about zero 
lift. Therefore, an algebraic expression of the drag 
coefficient in terms of C,? is not valid. This is partic- 
ularly true of polars of airplanes with deflected flaps. 
In order to have a general expression applicable to any 
polar, the effective induced drag coefficient is expressed 
as (Cz — Crmin-)?/wA.R., in the method presented here. 
This furnishes an easy means of approximating any 
polar with a negligible error. 

The assumption of constant parasite drag coefficient 
and the use of the airplane efficiency factor used in 
references 1, 2, and 3 are incorporated in the method 
shown here. However, the constant drag consists of 
the minimum drag plus the minimum induced drag. 
Because of the fact that the minimum drag does not 
necessarily occur at zero lift, the minimum induced 
drag is the induced drag resulting from the lift at 
minimum drag. When (Cy — Cymin,.)? is expanded, a 
linear function of C, also appears in the expression for 
the drag coefficient. In the final performance equation 
this term is accounted for by use of the coefficient of 
asymmetry. 

The basic performance equation is developed by use 
of the loading parameters and used in refer- 
ences 1, 2, and 3. The velocity at maximum lift-drag 
ratio is used herein as the basic velocity. It is expressed 
in terms of the loading parameters \, and d,. Hence, 


2 
WW 
ro) 
= 
DRAG COEFFICIENT, Cop 
Fic. 1. Typical airplane polar. 


DRAG COEFFICIENT, Cp 


Fic. 2. Drag coefficient of a PL airplane as a function of 
(Cr — Cr. for Cos = Oand Cr Lmin, = 0-4. 


it is independent of the power and propeller efficiency, 
Velocity is expressed in terms of the basic velocity by 
use of the nondimensional velocity ratio. The rate of 
climb is found at a predetermined velocity. This 
eliminates the necessity of assuming a _ propeller 
efficiency. However, determination of level flight 
velocities requires an assumed propeller efficiency, and 
two or three trials may be necessary in order to obtain 
an exact speed. 

Because of the various types of power plants, the 
complex nature of propeller efficiencies, and the use of 
exhaust jet thrust on modern airplanes, no attempt is 
made to obtain an analytical expression for power in 
terms of velocity. 

It is believed that a rapid performance evaluation 
with a high degree of accuracy may be achieved by the 
use of this method. 


DISCUSSION OF POLAR 


According to the Prandtl wing theory, the induced 
drag of an elliptical wing is proportional to the square 
of the lift coefficient. However, because of the lift and 
drag contributed by the tail plane and the variation of 
profile drag with angle of attack, the effective induced 
drag, including the variation of profile drag with angle 
of attack, is more closely approximated as a linear 
function of (C, — Crmin.)®. Essentially, the polar is a 
parabola that is asymmetrical with respect to the lift 
coefficient. A typical polar is shown in Fig. 1. It can 
be seen that the minimum drag does not occur at zero 
lift. Therefore, the approximating parabola is not 
symmetrical about zero lift. The drag coefficient is 
expressed as: 


Co = Comin. + — (1) 
By expanding (C, — Crm)? and substituting Cp. = 


Crmin./ TA-Ree and Coy = Coste. Eq. (1) re- 
duces to 


Cp = Coy + (C,?/mA.R..) 2CpaC, (2) 
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Hence, the polar for any airplane may be defined by 
Comin. 200 Crmin, and the airplane efficiency factor e if 
the geometry of the airplane is known. The value of e 
js easily determined by plotting the drag coefficient 
obtained from test data against (Cz — Cymin.)®. Crmin. 
and Cpmin. are easily determined from observation of the 
polar. Crmin. is the value of lift coefficient at minimum 
drag coefficient. A plot of drag coefficient versus 
(C, — Crmin.)” for a typical airplane with flaps deflected 
is shown in Fig. 2. The flapped configuration is used 
because it best illustrates this method. The drag 
coefficient is plotted against (Cz — Cymin.)? for Comin. = 
0 and Crmin. = 0.4. It is readily seen that a close 
approximation of the test polar is obtained by using 
Cimin. = 0.4 for this particular airplane. Observation 
of the curve for Crmin. = 0 shows that the analytical 
expression of the drag coefficient in terms of C,? does 
not approximate the polar so accurately as the curve for 
= 0.4. 


DERIVATION 


The determination of performance, regardless of the 
method, originates with the energy equation parallel to 
the flight path. 


WV sin @ = W(dh/dt) = TV — DV (3) 


Consider the drag term first. From Eq. (2) 
Co = (f/s) + (4) 
Therefore, 
207 
D 2CpaW (5) 
fpV? 2Ww? 
= —2 
DV 9 + CoaWV (6) 
The term TV is obviously (b.hp.)(nA). Hence, the 
energy equation becomes 
‘ 
W(dh/dt) = (b.hp.)n xpVb% +2Cp,.WV (7) 
and the expression for rate of climb is 
(b.hp.)n 2wW 
dh/dt = ———— —— 


In Eq. (8) p is replaced by poo and V is replaced by 
Vo. 


dh fooV? 2wW 

dt W 2w apo V + 
(9) 

Basic Speed 


The speed at the maximum value of L/D is now de- 
termined. D/L is expressed in terms of C;. 


D Cp Coy Cz 


(10) 


D/L is differentiated with respect to C, and equated to 
zero. The lift coefficient at minimum D/L is deter- 
mined by solving for C,. This obviously equals the lift 
coefficient at maximum L/D. Therefore, at the maxi- 
mum value of L/D, 


C, = VxA.R.Cp, (11) 


Hence, the indicated velocity at maximum L/D is 


Ve = V2W/pSC; = V2W/poV. rb’ef (12) 
Substituting the loading parameters in Eq. (12), 
j 4 
Ve = (+) 1,1/4],1/4 (13) 
Hence, Vz = YT’. In order to express any velocity in 
terms of T’, the parameter R, is introduced. 
Vs = 


If the loading parameters dy, d,, and d, are substi- 
tuted into Eq. (9), the equation for the rate of climb 
becomes 


3 
Since any velocity may be represented by the basic 
speed, Vz = d,'“d, “‘ and a velocity ratio R,, 


dt Ap ( R, 

(15) 


From the definition of C, it can be seen that 


Therefore, 

dh 1 2C, 
t dt 1 vi (R, + R, a 


1 =l]- (16) 


Level Flight Velocity 

The required condition for level flight velocity is zero 
rate of cimb. The relationship between the basic per- 
formance parameter and tle velocity ratio is obtained 
by setting dh/dt = 0 in Eq. (16). Therefore, in level 
flight, the relationship between the basic performance 
parameter and the velocity ratio is 


= 1/[R® + (1/R,) — (17) 


Time to Climb 


The time required to climb to any altitude with con- 
stant power is determined by integrating the reciprocal 
of the rate of climb with respect to altitude as shown 
below. It is shown in Eq. (16) that the product of the 
rate of climb times the power loading is (1 — y,’R,). 
Therefore, the reciprocal of the rate of climb is 


1/(dh/dt) = d,/(1 — (18) 
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Fic. 3. Chart for determination of basic velocity and rate of 
climb. 


The integral of the rate of climb with respect to altitude 
is 


T. = [1/(1 — (19) 


vy,’ is replaced by y’/+/c in Eq. (19). Hence, the 
value 7./, may be expressed as 


m 1— m Vo 
Effect of Gross Weight on Climb 
Since \, is a linear function of weight and y is a func- 
tion of W, the general energy equation may be written 


in terms of an assumed weight and the weight ratio R,, as 
shown helow: 


R,.(dh/dt) 2y,/RyRw (21) 
SIGNIFICANCE OF PARAMETERS 


Basic Performance Parameter 
The dimensionless parameter y;’ is indicative of the 
airplane performance. However, the parameter yp,’ is 
the true index of climb performance at a specified ve- 
locity. A perfect airplane is represented when the basic 
performance parameter is zero. Although this value 
may only be achieved by infinite power, infinite span, 
or zero parasite drag, it shows that the airplane climb 
performance improves as ~g becomes ‘smaller. The 
value \,(dh/dt) = 0 and the absolute ceiling of an 
airplane occur when yp,’ equals 0.5. An airplane with 
ve’ equals 0.5 at sea level will not fly. If the perform- 
ance parameter is expressed in terms of feet, pounds, and 
horsepower, ¥z will equal 22.37 when pp’ equals 0.5. 
Concerning the level flight speed, the value y,’ is the 
index of performance as Ry is a function of the inde- 


JOURNAL OF THE AERONAUTICAL SCIENCES—APRIL, 1947 


pendent variable R,. The velocity in level flight in. 
creases as the basic performance parameter becomes 
smaller. A theoretically perfect airplane with in. 
finite velocity is represented by the value of y being 
zero. 

The basic performance parameter y has a definite re. 
lationship to the fundamental performance parameter A 
used in references 1, 2, and 3. From the definitions of 
y and A it is readily seen that y = A’*. There are 
several reasons for the use of ¥ in preference to A in the 
method presented here. The derivation of the formulas 
evolves into a final expression that is a linear function of 
y. The relationship of R, to R, is proportional to the 
relationship of speed to power. Hence, observation of 
Fig. 4 permits visual determination of the position on 
the power required curve and R, has a definite physical 
significance, 

Observation of Fig. 3 shows that /,(dh/dt) is a linear 
function of y. Therefore, the value of dh/dt may be 
readily determined from ~r without the use of charts if 
desired. 


Basic Speed Parameter 

The basic speed parameter T furnishes a useful tool in 
performance determination. It is numerically equal to 
the fourth power of the velocity for maximum L/D. In 
terms of the known airplane parameters, it is A,A,. It 
is the only parameter used herein that is dimensional. 
It has the dimensions of (length/time)*. The actual 
units depend, of course, on what units the performance 
is expressed in. 


The Parameters R, and Ry 

The relationship between R, and R, is obviously the 
relationship between speed and power. The velocity 
ratio, R, is the ratio between any desired speed and the 
speed at which maximum L/D occurs. R, is the ratio 
between the basic performance parameter based on any 
desired speed and the basic performance parameter 
based on the speed for maximum L/D. Inasmuch as 
plot of R, vs. R, is a nondimensional speed power curve, 
a desirable climb speed may be determined from 
observation of the curve. This may be seen from 
observation of Fig. 4. 


The Coefficient of Asymmetry 


The coefficient of asymmetry is indicative of the 
symmetry of the polar. At C, equals zero, the polaris 
symmetrical about zero lift. For any asymmetrical 
polar C, has a finite value. This coefficient is incor- 
porated in the parameter R, as shown in Eq. (16). It 
should be noted that the asymmetry of the polar is 
accounted for in two ways. The constant term Cpimin. 
is added to the minimum drag coefficient to give a fie- 
titious equivalent minimum drag. The linear function 
of lift coefficient is accounted for by use of the non- 
dimensional coefficient C,. 
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Fic.4. Performance parameter ratio as a function of climb speed 
ratio and Ca. 
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Fic. 5. Level flight velocity ratio as a function of y and Cuz. 
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Fic.6. Chart for determining time required for a constant power 
climb. 


DISCUSSION OF CHARTS 


Charts for Determining Rate of Climb 


In order to determine rate of climb, it is necessary 
to use the chart showing the performance parameter 
ratio as a function of the climb-speed ratio shown in 
Fig. 4. This chart consists of a family of curves for 
various values of C,. It should be noted that /,(dh/dt) 
may be easily calculated for any values of y, C,, and R,. 
Inasmuch as no charts are required for the determina- 
tion of the climb-speed ratio and the determination of 
the rate of climb, the charts shown in Figs. 2 and 3 are 
not essential for determination of rate of climb. Alge- 
braic solution of Eq. (16) is advocated if an exact value 
of rate of climb is desired. However, the charts furnish 
a reasonably accurate method to determine rate of climb. 


Charts for Determining Velocity 


Because of the relationship between R, and R, being 
a fourth degree equation, the complexity of an algebraic 
solution is prohibitive. Therefore, it is desirable to 
determine velocity by use of the charts. By assuming 
values of R,, a family of curves with various values of 
C,, may be plotted as shown in Fig. 5. This furnishes a 
convenient manner of obtaining level flight velocity. 
An interesting characteristic of this chart is that it may 
be used in reverse. If a specified speed is desired, the 
chart may be used to determine the required power by 
determining from the chart the value of ¥ required for 
the specified value of R, and solving for the required 
thrust horsepower. 


Charts for Determination of Time to Climb 


It is difficult to establish a general chart to determine 
the time required to climb because of the infinite possi- 
bilities of application of power. However, with the 
advent of torquemeters in modern airplanes, climbs at 
constant power are now of practical significance. Fig. 6 
shows 7,/l, plotted vs. y for various altitudes. This 
chart furnishes a convenient means of obtaining time 
to climb to any specified altitude at a constant power. 


Effect of Gross Weight on Climb 


Heretofore, the performance has been considered at a 
constant assumed weight. In determining rate of 
climb, it is often necessary to know the weight that will 
allow a specified climb for a given airplane configura- 
tion. For transport operation, this is particularly im- 
portant, since the Civil Aeronautics Administration 
specifies certain rates of climb for various airplane con- 
figurations. For military operation, it is important to 
know at what weight the absolute ceiling or the service 
ceiling will occur at a specified altitude. Fig. 7 shows 
the weight ratio plotted versus the performance par- 
ameter for various values of /,(dh/dt). From this chart 
it is easy to determine the change in weight from an 
initially assumed weight required to attain a specified 
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Fic. 7. Weight ratio for a specified climb as a function of basic 
performance parameter. 


rate of climb. By using the curve for /,(dh/dt) equals 
zero, the change in weight required for absolute ceiling 
is obtained. 


APPLICATION OF METHOD 


In applying this method to actual use it is necessary 
to determine T, ¥;,, and C, for a basic weight at sea 
level. The parameter ¥;, may be determined at any 
weight by multiplying by the factor (W/Whasic)’’, and 
it may be determined at any altitude by multiplying 
by the factor ~/ p/p. In order to determine y, at some 
power other than the assumed power, y, is multiplied by 
the ratio of the new power to the original assumed 
power. The parameter T is independent of both power 
and altitude and is directly proportional to W?. The 
parameter C, is a function of the polar and, once 
established, does not change unless the airplane con- 
figuration is altered. In actual application, the engi- 
neering units are used and the necessary conversion 
factors are used in determining the formulas and charts. 
The following units are used to determine the parame- 
ters in the charts shown herein: 


1, = Ibs. per sq.ft. 
l, = Ibs. per sq.ft. 
l, = Ibs. per t-hp. 


Determination of Rate of Climb 

The rate of climb may be determined by use of 
formulas or by use of charts. The procedure with 
charts is outlined below: 


(1) Determine the value of T, y,, and C, for a specified 
configuration and altitude. 

(2) Obtain the basic speed from Fig. 3. 

(3) Obtain the value of R, from Fig. 4 as a function 
of R,and C,. If the maximum rate of climb is desired, 
a value of R, equals 0.9 will closely approximate this 
for a conventional constant speed propeller and a re- 
ciprocating engine. If the exact maximum value is 
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required, it is necessary to choose three values of R, and 
plot /,(dh/dt) versus R, to obtain a maximum value. 
(4) Correct y, for desired velocity by multiplying by 


(5) Obtain a value of /,(dh/dt) from Fig. 3 and divide 
by /, to obtain the rate of climb. 


Level Flight Velocity 


The procedure used to obtain level flight speed is 
shown below. 

(1) Assume a propeller efficiency and obtain the 
value of T, ¥,, and C,. 

(2) Obtain V; as a function of T from Fig. 3. 

(3) Determine R, from Fig. 5. Velocity is R,(V,) X 


po/). 


If the assumed propeller efficiency is incorrect for the 
velocity found in step 3, determine a new propeller 
efficiency and a new y, by multiplying original y, by 
n/n and repeat step 3. One or two repetitions usually 
provide an accurate answer. 


Time to Climb 


The procedure for determining the time to climb to 
any altitude at constant power is outlined below: 


(1) Obtain the value of y, at sea level. 
(2) Obtain from Fig. 4 the value R, and multiply y, 


by R, to obtain yp. 
(3) Determine T,/l, from Fig. 6 and obtain the time 
to climb by multiplying T,//, by /,. 


Effect of Weight on Climb 

The procedure for determining the change in weight 
necessary to obtain a specified rate of climb is shown 
below: 

(1) Determine y, and R, and /,(dh/dt). 

(2) From Fig. 7 obtain directly the weight ratio re- 
quired to obtain the rate of climb specified in step 1 and 
multiply the assumed weight by the weight ratio to 
determine the weight that permits the specified rate of 
climb. 


NONDIMENSIONAL SPEED POWER CURVE 


From observation of Eq. (16), it can be seen that 


R, = 1/2y,’ during level flight. The relationship be- 
tween y,’ and Pj is: 
1 P 
=—“=R, (22) 
2y, As Asp 


Hence, the relationship between R, and R,,_, is: 
Row = + 
A plot of R,,, vs. R, is shown on Fig. 8. This curve has 


a unigue application, since it is a nondimensional speed 
power curve and furnishes a means of determining 4 
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graphical or trial-and-error solution for level flight 
velocity and rate of climb. The use of this curve per- 
mits observation of the position on the power curve 
without the laborious process of computing power re- 
quired versus velocity for various altitudes and weights. 

Since R,,, is equal to R,, Fig. 4 differs from Fig. 8 
only in that the former is plotted only over the range of 
climb speeds and the latter is plotted over the range of 
speeds from 80 per cent of speed for maximum L/D to 
high level flight speeds. 

From the energy equations shown in Eqs. (7) and (8), 
it can be seen that 


AP 


WAW/W) 


dh/dt = AP/W = 


From the definition of Pi, dh/dt may be expressed as 
shown below: 


dh _ 
dt W, Vo 
Hence, the value of Cj, is: 


AR,,, "A 
W, 


(25) 


Ciw (26) 


From the above relationships, it can be seen that the 
rate of climb is readily determined from Fig. 8 by plot- 
ting the values of R,,, vs. R, for any specified values 
power and altitude. 

The procedure for determination of performance is 
shown below. A standard weight, W,, is first assumed. 
Preliminary to determining performance for specific 
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Fic. 8. Nondimensional power required curve. 


altitudes, weight, and powers, it is necessary to deter- 
mine the following values: 
(1) Determine T and V;,,, for a standard weight, W,. 
(2) Obtain the value of C, and },,“* for a standard 
weight, WV,, 


= 0.0447(W,/b Ve)" 


The following shows the procedure for determining 
level flight speed for any condition: 


(1) Assume a velocity, V. 

(2) Determine R, = V/Vim V(W/W,)o. 

(3) Determine 7 as a function of power, engine speed, 
velocity, and propeller diameter. 

(4) Determine Rp», = (b.hp.)n 
and plot vs. R, on Fig. 8. 

(5) Repeat above procedure if the point obtained in 
step 4 does not fall on the curve for the value of C, de- 
termined above. Two repetitions of the above pro- 
cedure furnish enough points to plot a curve of Rp». vs. 
R, on Fig. 8. The value of R, for level flight is deter- 
mined from the intersection of the calculated curve and 
the curve for the correct value of C, on Fig. 8. 

(6) Determine V = R,Vz,, Vv ( W/W,)/ for level 
flight. 


The procedure for determination of maximum rate of 
climb for any condition is shown below: 


(1) Assume several values of velocity and compute 
n as a function of power, engine speed, velocity, and pro- 
peller diameter. 

(2) Determine the values of R, for the values of 
velocity shown above. 

(3) Compute R,,,. for the values of 7 and V in step 1 
and plot vs. Ro. 

(4) Measure the maximum value of AR, up 10T the 
previously determined value of C, and compute C,, = 
33,000. 

(5) Compute dh/dt = Civ V W/W,/ Vo ft. per min. 

It should be noted that the above procedure must be 
completed for any weight. The value of AR, and Cj» 
in steps 4 and 5 are for a specific weight, and dh/dt may 
not be determined for other weights without determin- 
ing a new value of P,,,. 


CONCLUSIONS 


By mathematically expressing the airplane polar 
about an axis of symmetry that coincides with the 
minimum drag coefficient, taken from wind tunnel or 
flight data, airplane performance may be determined 
analytically with negligible error. 

The method of performance presented above furnishes 
a convenient manner to determine rate of climb, level 
flight speed, and time to climb. The rate of climb is 
obtained for a predetermined velocity, eliminating the 
necessity of assuming a speed and a propeller efficiency. 
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Changes in performance due to changes in weight or 
power are easily determined. The change in weight 
from an assumed weight to obtain a specified rate of 
climb may be determined directly from a chart. 

A performance parameter is presented which is in- 
dicative of the performance of an airplane and furnishes 
a qualitative comparison of airplanes without a detailed 
analysis. 

Although an assumed propeller efficiency is necessary 
to determine level flight velocities, a trial-and-error pro- 
cedure is simple because of the basic performance 
parameter being directly proportional to the propel- 
ler efficiency. 


A nondimensional power required curve is presented 
for graphical performance determination without the 
laborious process of computing power required curves 
for various weights and altitudes. 
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Letters to 


Dear Sir:. 
Starting with the general stress-strain curve equation given by 


Ramberg and Osgood! 
s 
E +K (1) 


and using as parameters the basic modulus of the material (Z) and 
two points not on the straight line portion of the actual stress- 
strain curve, and F2) rather than the Ramberg-Osgood 
stresses Foz and Foss; (where the subscript denotes stress at a 
secant modulus to tangent modulus ratio of 0.70 and 0.85, respec- 
tively), the following equation is obtained for the stress-strain 
curve: 


e = (F/E) + [a — (Fi/E)](F/Fi)" (2) 
where ¢ is the desired strain at stress F, and m is obtained from 


log F./F; 


where ¢; = strain at stress Fj, = strainat stress F2,,and F, > F,. 
Continuing in this fashion, the equation of the (stress) versus 
(stress/tangent modulus) ratio curve becomes: 


F F 
(4) 


Now, since 
e = F/E, (5) 


where E, = secant modulus. 
The substitution of Eq. (5) in Eq. (2) results in: 


F F Fi\f/ F \" 
i + (« (6) 
‘ Comparison of Eqs. (4) and (6) shows that the (F/E,) versus 
(F) curve may be obtained directly from the material stress- 
strain curve by the simple method of multiplying the offsets of 
the original stress-strain curve from the basic straight modulus 
line by the coefficient ” and plotting through the points thus ob- 


tained. 
Cozzone and Melcon in their paper? point out that equations of 


the form 
F, = w*E,/(L'/p)? (7) 
can be rearranged so that 


E./F. = (L'/p)*/x? (8) 


the Editor 


If we now invert Eq. (8), we obtain 
F/Ex = w*/(L'/p)* = cx*/(L/p)? (9) 


where c = column end fixity coefficient. 
If we let 


EAR (10) 
then the term A takes on the significance of a strain. 

Now, since every value of A corresponds to a particular value 
of cr*/(L/p)?, it is possible to determine the column curve of the 
material from the stress-strain curve of that material by means of 
the following steps: 

(1) From a conventional stress-strain curve of the material 
determine n from Eq. (3) by measuring the strains ¢ and « at 
stresses F, and F2, respectively, and the basic slope of the stress- 
strain curve. Note that F; must be greater than F;. 

(2) Multiply the offsets of the stress-strain curve from the 
basic modulus line by the coefficient m determined in step (1) 
above and plot the F versus F/E; curve on the same sheet (if 
desired) with the original stress-strain curve. 

(3) From the relationships expressed in Eq. (10) compute the 
value of A for each value of L’/p. 

(4) Enter the F versus F/E; curve plotted per step (2) above 
at the value of A corresponding to each value of L’/p and proceed 
up vertically to the F versus F/E; curve, thence to the left and 
obtain the values of F, corresponding to the value of A = cz*/- 
(L/p)?. 

(5) From the values of F, and L’/p thus obtained the column 
curve of the material may easily be obtained. 

This method may be extended to other allowable stresses in- 
volving the tangent modulus £; or the secant modulus E£, after 
the method of Gerard. 

Fhis method of determining column curves is also applicable 
to the Ramberg-Osgood equations,! since with their general nota- 


tion 
Ss — m\fs \" 
- 11) 


after proper differentiation 


s\f1 — m\f s \" 
12) 


where 7 is evaluated by Fig. 10 of reference 1. 
The relationship between Eqs. (11) and (12) is precisely the 
same as that existing between Eqs. (2) and (4). The only differ- 
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ence exists in terminology and in the fact that with Eq. (2) any 
two points on the stress-strain curve not on the straight-line por- 
tion thereof may be used to determine the coefficient , rather 
than stresses at specific values of the E,/E ratio. 

Other relationships which result from operating upon Eq. (2) 


are: 


EF," (13) 
+ Ba — 

(14) 
F\" 

E, = EF,"/[F." + F"- (Ea — F,)} (15) 
E, 1 
(16) 


If the proportional limit stress Fp, at a particular value of 
permanent set or deviation from the straight-line portion of the 


stress-strain curve is desired, then x 
For, = E( dpa.) /(Ea — Fi) (17) 


where Ap,, is the desired limit for deviation from the straight line 
stress-strain curve. 
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LEE Gray 
Stress Project Engineer 
Wichita Division 
Boeing Airplane Company 


Dear Sir: 

The paper by A. H. Flax appearing in the January, 1947, 
issue of the JoURNAL presented an excellent summary of blade 
stress analysis procedures but seems to have neglected several 
significant factors that might well have been considered. 

(1) The writer has shown’? that the basic equation govern- 
ing blade bending is reducible to a linear second order equation 
with variable coefficients, the equation being as follows:? 


x x 
+ Ele” = SS f(x)dx — ByS(sin 8 + ¢ cos 8) S, Auxdx 
where 


EI = product of modulus of elasticity times moment of 
inertia (variable) 

product of cross-sectional area times unit weight 
(variable) 

S = rotor semispan 

Bo = w cos 3/2; w = rotational velocity 

8 = precone angle or angle between undeflected blade axis 

and plane of rotation 


Il 


= d{/dt; = S(dy/dt) = de/dx; = = 

£ = undeflected blade axis 

= axis normal to 


x 
f(x) = arbitrary load distribution which need mot be ana- 
lytically definable 


The boundary condition being ¢(0) = ¢’(1) = 0 for blades fixed 
at the root or continuous across the rotor mast, ¢’(0) = ¢’(1) = 0 
for hinged blades without dampeners, and ¢g’(0) = SM;/EI,, 
¢’(1) = 0 for hinged blades with dampeners, where 14, = hinge 
moment and J, = root moment of inertia. 

The two additional available boundary conditions—namely, 
¢(0) = 0 and shear at tip = O—have already been used in the 
definite integrals used to define the shear and deflection. Note 


x 
that ¢ = deflection = SJ, dx (since ¢ is defined by the basic 
equation this latter integral is not used to determine ¢). 

It should be noted that the above equation is readily modified 
to include effects of items such as an offset hinge axis. 

(2) Depending on what set of boundary conditions are 
involved, the above equation may or may not be soluble by the 
various methods outlined in Mr. Flax’s paper. Experience tends 
to indicate that the hinged rotor case is by far the easiest to 
handle, and perhaps that is why considerable literature is avail- 
able on the subject. When dealing with rotors with fixed ends 
or rotors continuous across the mast (such as found on Bell, 
Bendix, Hiller, and other machines), it appears that only one of 
the methods mentioned in Mr. Flax’s paper will work. This is 
the Stuart method.* In addition to the Stuart method, the other 
available means for solving the equations of blade bending is 
the one utilizing finite differences as outlined in the writer’s 
papers.!? 

(3) The particular advantage of the finite difference method 
is that under no circumstances are any guesses necessary as to 
what the root moments should be, etc. The method is direct 
and leads at once to the complete determination of blade loads 
in one operation. Furthermore, like the Stuart method, it is 
applicable to any type of rotor blade subjected to any type of 
loading, since a semigraphical solution is readily obtained. 

It is regrettable that all methods of blade analysis today 
present a rather terrifying appearance to all but those familiar 
with mathematical procedures. Invariably one gets a first 
impression of tremendous complication and of a tremendous 
amount of work involved. Perhaps more time should be spent 
in discussing physical principles in an attempt to remove the 
terrifying aspects of Greek letters and mathematical notations. 
Then, too, more time should be spent in describing the actual 
operations involved in actually obtaining numerical results. 
In this regard the writer has been asked a number of times how, 
in a reasonable length of time, is it possible to solve the five, 
ten, or 20 simultaneous equations involved in the finite difference 
solution for blade bending? When one is familiar with the 
equations involved it is clear that, because of the nature of the 
equations, even 20 of them can be solved in less than 8 man-hours. 
When only five are involved, work of solution is almost insignifi- 
cant compared with the effort expended in computing the blade 
section properties. This rapid and simple solution is, of course, 
not a characteristic of simultaneous equations in general, and, 
unless the particular equations are examined in detail, it is not 
at all obvious that little work is involved in their solution. The 
reason for simple solution in the case under discussion is simply 
that none of the equations contains more than three of the 
unknowns and that these unknowns are grouped ascendingly. 
In addition, the boundary conditions remove an unknown from 
the first and last equations so that a system is obtained as follows: 

If gn = g at any point, then, if a fixed blade be considered, 
g = 0 and the following typical system results: 


g 
Rsgi — €3 
= Risgi 
= — €19 


with Ciga + Cogis — Creo = 0, where Ci, C2, and C; are constants 
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But ¢i9 = ¢a from the boundary conditions, so that 


(Ci + Cr)as — 


where the & and e terms are constants depending on blade 
properties, loading and so forth. Most engineers have solved 
systems of equations almost identical to the above in connection 
with bolt pattern analysis, yet experience has shown that unless 
the similarity is pointed out, it will not be recognized. 


SUMMARY 


(1) In addition to the blade analysis methods discussed in 
Mr. Flax’s paper, the method of finite differences is available as a 
general rapid means to solve the blade bending equations. Finite 
differences have the advantage that they may be applied to any 
known type of rotor system and do not demand that blade shape 
or loading be definable analytically. 

(2) The basic equation governing rotor bending moments is 
a second order differential equation with the slope rather than the 
deflection as the dependent variable. 

(3) Only the Stuart method in addition to finite difference 
methods has been shown to yield solutions for the structurally 
continuous or the fixed end blade system. 


(4) Although this is not intended to apply particularly ty 
Mr. Flax’s paper, since the writer himself is far from blameless, 
it is felt that authors of technical papers in general could ensure 
far greater utility of their products by endeavoring to show 
clearly the advantages of any new procedure, not only from the 
point of view of brilliant technical discovery but also from ay 
application standpoint. Since the most popular procedures 
are the ones that we understand best and which we can apply 
easily, it would seem that more emphasis should be placed op 
giving the reader a good physical picture of the method, as wel] 
as an understanding of how the numerical drudgery can be 
minimized. 
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DANIEL O. DoMMASCH 
Engineering Coordinator 
Bell Aircraft Corporation 
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Shear Distribution 


Due to Twist in a 
Cylindrical Fuselage with a Cutout 


M. M. GOLOLOBOV* 
Letov, Prague 


SUMMARY 


A general method is given by which the shear distribution due 
to twist in a cylindrical fuselage with a cutout can be obtained. 
In addition, a simple example is solved to illustrate the shear dis- 
tribution. 


INTRODUCTION 


VERY FUSELAGE has several cutouts, the largest, of 
course, being those for cabin doors. These cut- 
outs cause considerable changes in the distribution of 
shear stresses due to torsion and additional axial forces 
in the stringers. This paper presents a method of 
solving the stated problem. 

In order to obtain a simple solution of this statically 
indeterminate problem, the following assumptions are 
made:f (a) the stringers plus effective width of skin 
carry only axial forces; and (b) the skin between the 
stringers is subjected to shear stresses only. 


r GENERAL METHOD 


Let us consider a cylindrical fuselage of circular cross 
section with a rectangular cutout (see Fig. 1). On both 
sides of the cutout are reinforcing rings. In the follow- 
ing it is assumed that the changes in the shear distribu- 
tion caused by the cutout are confined to the length 
shown in Fig. 1, i.e., that they do not reach farther than 
the two adjacent bays of the fuselage. 

First of all, let us determine the shear distribution 
in the section through the middle of the cutout (section 
D-D in Fig. 1). This can be done in the following 
manner: a shear g, per unit length is applied along the 
whole perimeter of the circle (see Fig. 2); then a force 
@> (equal to the resultant of a shear —q acting along 
the cutaway perimeter) is applied to cancel the shear 
along the cutout, distance h being h = 25S,/b. As the 
resultant shear distribution must give only a torque, it 
is further necessary to apply a force of the same magni- 
tude but in the opposite direction at the shear center 
C of the section with the cutout. 

The shear g, is found from the equation T = 2Sq, — 
@bc so that 


Received June 11, 1946. 

* Chief Stress Analyst. 

} These assumptions are usually made when solving statically 
indeterminate structures. The real stress distribution differs 
from this idealized one in that the skin is subjected not only to 
shear but also to direct stresses. 
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Q@ = T/(2S — be) (1) 


where S = 2D?/4. 

The shear distribution in the middle section is given 
by the algebraic sum of the uniform shear q and the 
shear distribution due to the shear force qb acting at 
the shear center C, it is therefore necessary to deter- 
mine this shear distribution. 

The shape of the middle section is shown in Fig. 3. 
The two stringers at the edge of the cutout are stronger 
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than the rest; let us denote their cross-sectional area 
A;. The area of the cross section of all other members 
effective in bending stringers pius effective width of 
skin, may be replaced by a fictive skin of thickness s,. 
It is a well-known fact that the shear at any point E of 
the cross section due to a shear force F is given by 
q = FQ/I, where Q is the static moment of the cross- 
sectional area effective in bending between points D 
and E about the neutral axis and J is the moment of 
inertia of the whole cross section effective in bending. 
The shear per unit length at point E defined by angle ¢ 
(see Fig. 3) due to a shear force F = 1 lb. applied at the 
shear center C is 


Ad np 
= or +7 rdgsr sin ¢ 
or 
r’s, 
qe = oF + va (cos a — cos ¢) (2) 


From Eq. (1) it is evident that we must know the 
location of the shear center to be able to determine the 
shear g,. The sum of the moments of the shears q, 
about the center of the section O is given by: 


M = Sq,rder = 


an f o (cos a — cos g)dg = If (3) 


where f is the distance of the shear center C from the 
center O of the cross section (see Fig. 2). The location 
of the shear center may be found without evaluating 
the integrals by the following reasoning: The first 
term of Eq. (3) represents the moment of a constant 
shear A,b/2/ the resultant of which is F; = A,b?/27 and 


Fie. 5. 


acts at a distance k’ = 2S,/b from the edges of the cut- 
out (see Fig. 4). The remaining shear force is F,; = 
1 — A,b?/2I = I,/I where J, is the moment of inertia of 
the cross section without the edge stiffeners. The sec- 
ond term of Eq. (3) may therefore be represented as: 


A,b? 
2r 7 (cos a — cos (1 y 


where / is the distance of the shear center of the circular 
are from point O. 
This distance is given by the expression 


= 2r[(sin — 8 cos B)/(8 — 1/2 sin 28)] (4) 
The distance f is then 
1f = (Ayb?/2Dk + [1 — 
or 
f = 1+ (Aib?/21)(k — (5) 


Now we can determine the shear gq, from Eq. (1). 
The resultant shear will be given by the expression (see 
also Eq. (2)): 


qo — | 


r’s, | 
q2 


or t+ 7 — 


or 


q2 = T cos @ cos ¢ 


Substituting for g, from Eq. (1) we get, finally: 


A,b? br?s, 


The resultant shear distribution is shown in Fig. 5a. 
As we have replaced the stringers by a skin of thick- 
ness s, we have got a smooth curve. The real sheaf 
distribution differs somewhat from this idealized one: 
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between two neighboring stringers the shear remains 
constant and there is a sudden change of shear at every 
stringer. This shear distribution is shown in Fig. 5b. 
The difference of shears to the left and to the right of 
every stringer sets up an axial force init. The diagram 
of the axial force in a stringer is shown in Fig. 6. The 
axial force at the middle section is zero; it increases 
linearly up to its maximum value at the section at the 
edge of the cutout, and then diminishes linearly to zero. 

The shear distribution in the part of the fuselage ad- 
jacent to the cutout can be determined in the following 
manner. The difference Aq, of shears to the left and 
to the right of a certain stringer in bay I (see Fig. 1) 
can be found from the difference Ag, of shears to the 
left and to the right of the same stringer in bay II from 
the equation 


Po = Aqe = Aq(d/2) 
so that 
Aq = (d/2e) 


It follows, therefore, that the shear distribution in 
bay I is obtained by multiplying the shear distribution 
in bay II (with the cutout) by d/2e and adding a con- 
stant shear Ag of such magnitude that the resultant 
torque at this section becomes equal to JT. This con- 
stant shear Ag is given by: Ag = [1 — (d/2e)](T/2S). 

The resultant shear distribution in bay I is therefore 


d 
(7) 


or 


E Ab? br’s, 


= 3025 — be (cos a=- cos + 


(1 3) 

2e/ 2S 
The resultant shear distribution in bay I is shown in 
Fig. 7 for the case that d/2e < 1. 


I 


NUMERICAL EXAMPLE 


As an illustration of the method let us determine the 
shear distribution due to torsion in the part of the fuse- 
lage shown in Fig. 8. The fuselage has a diameter of 
98.5 in.; the rectangular cutout is 52.75 in. by 27 in. 
The edge stiffeners have a cross-sectional area A; = 
0.93 sq.in.; stringers (inclusive of effective width of 
skin), A, = 0.194 sq.in. Stringer spacing is 6.18 so 
that s, becomes s, = 0.194/6.18 = 0.0314. The 
moment of inertia of the fuselage cross section (with 
the cutout) about the axis X-X is J = 12,380 in.‘ 

First of all we must find the location of the shear cen- 
ter of the section with the cutout (see Figs. 2 and 4). 


S = rD?*/4 = 7,620 sq.in. 
S, = 270 sq.in. 


Fic. 6. 
Fic. 7. 
X 
= 
D=98> 
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270 
= = —— 93 in. 
h 2S,/b = 2X 52.75 10.23 in 


Sx = S— S, = 7,350 sq.in. 
Distance k’ (see Fig. 4) 


2 X 7,350 
= = = 279 in. 
k 2S,./b 59.75 79 in 
k = k’ — (98.5/2 — 7.67) 
k = 237.4 in. 
6B = 180 — @ = 147.5° 
sin 8B — Bcos B 
= 2 = 98.5 
: B — */ssin 28 = 
0.536 + 2.574 XK 0.844 . 
= 88.2 in. 
2.574 + 1/2(0.906) 
= 88.2 in. 
From Eq. (5) 


0.93 X 52.75? 
2 X 12,380 


(237.4 — 88.2) 


2 
14+ = 88.24 


 Cut- 
| 
= 
tia of 
Sec- 
as: 
(4) = | 
hick- 
sheaf 
= 


103.8 in. 

Distance c (see Fig. 2) 
= 103.8 + 41.54 + 10.23 
155.6 in. 


Now we can determine the shear g, (see Eq. (1)) 
T T 


% ~ 95—be 2X 7,620 — 52.75 X 155.6 
q = 0.0001422 X T Ibs. per in. 


The resultant shear g2 in the middle section is given 
by Eq. (6) 


| 


95 — be a1 ii 


0.0001422 71 — 0.1044 — 


(cos a — cos °| = 


52.75 X 49.25? X 0.0314 
12,380 


q@ = 0.0000885T + 0.0000461T cos ¢ 
(see Table 1). 


The shear distribution in the middle section is shown 
in Fig. 5. The shear has its maximum value at the 
edge of the cutout. If the fuselage had no cutout the 
shear would be: 


T/2 X 7,620 = 0.0000657 X T 


The maximum shear at the edge of the cutout is twice 
as great. 

Now we can determine the shear distribution in the 
adjacent parts of the fuselage (see Eq. (7)): 


(0.844 — cos 


d d\ T 
In our example d = 27 in.; e = 18in.; d/2e = 0.75; 
hence, 
g = 0.75 X go + 0.25 X 0.0000657T = 0.75q + 
0.0000164T 


(see Table 2). 

The shear distribution in the parts of the fuselage 
adjacent to the cutout is shown in Fig. 7. 

Having determined the shear distribution, it is now 
an easy matter to calculate the axial force in each 
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stringer from the difference of shears to both sides of ity 


TABLE 1 
¢° 
32.5 0.0001274 T 
60 0.0001115 T 


TABLE 2 


¢° 
0 
32 
32 
60 0.0001001 T 
90 0.0000828 T 
120 
150 
180 


0.0000655 T 
0.0000529 T 
0.0000482 T 


— 


e=/8" 


e= 18" _ 


The axial forces of all stringers at the middle section 
are equal to zero. The largest axial forces are, of 
course, in the edge stiffeners. The diagram of the 
axial force in one of the edge stiffeners is shown in Fig 
9. 


Py = 0.00012747(13.5) = 
(0.00011207 — 0.00001647)18 = 0.00172F 


CONCLUSIONS 


The method shown is general in nature and can be 
applied with small modifications to fuselages of every 
cross section. It is extremely suitable for fuselages 
having cross sections with two axes of symmetry. fj 
however, the fuselage has a cross section with only omé 
axis of symmetry the calculations become complicated 
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